# ---- CODE ----
source("Code/functionsV1.R")

# ---- osmacCode ----
adult <- read.table("Code/adult.data", sep = ",")
test <- read.table("Code/adult.test", sep = ",", skip = 1)
test$V15 <- gsub("\\.", "", test$V15)
adult <- rbind(adult, test)
adult1 <- subset(adult, select = c("V1", # age
                                   "V3", # fnlwgt
                                   "V5", # edu
                                   "V12", # loss
                                   "V13", # hours
                                   "V15", # income
                                   NULL))
adult1$V15 <- as.numeric(adult1$V15 == " >50K")
income.glm <- glm(V15 ~ ., data = adult1, family = "binomial")
summary(income.glm)


# ---- osmac00 ----
set.seed(123)
# ---- osmac01 ----
X <- cbind(1, as.matrix(adult1[, -dim(adult1)[2]]))
y <- adult1$V15
set.seed(123)
AdpOptSubLog(X, y, r0 = 500, r = 1000, optmethod = "A", data = adult1, 
             covariate = "V1 + V3 + V5 + V12 + V13")


# ---- osmacUW00 ----
set.seed(123)
# ---- osmacUW01 ----
AdpOptUWLog(X, y, r0 = 500, r = 1000, optmethod = "A", data = adult1, 
            covariate = "V1 + V3 + V5 + V12 + V13")


# ---- osmacPS00 ----
set.seed(123)
# ---- osmacPS01 ----
AdpOptPosLog(X, y, r0 = 500, r = 1000, optmethod = "A", data = adult1, 
             covariate = "V1 + V3 + V5 + V12 + V13")

# ---- poissonCode ----
hour <- read.csv("Code/hour.csv")
hour1 <- subset(hour, select = c("workingday",
                                 "temp",
                                 "hum",
                                 "windspeed",
                                 "cnt",
                                 NULL))
hour.glm <- glm(cnt ~ ., data = hour1, family = "quasipoisson")
summary(hour.glm)

# ---- poisson00 ----
set.seed(123)
# ---- poisson01 ----
y <- hour1$cnt
X <- cbind(1, as.matrix(hour1[, -dim(hour1)[2]]))
AdpOptSubPoi(X, y, r0 = 200, r = 500, optmethod = "A", 
             delta.quant = 0.05)


# ---- QuantileCode00 ----
casp <- read.csv("Code/CASP.csv")
casp <- casp[, -which(colnames(casp) == "F3")] # F3 = F2/F1
fit.full <- rq(RMSD ~ ., tau=0.75, data = casp, method="pfn")
summary(fit.full)

# ---- QuantileCode01 ----
set.seed(123)
# ---- Quantile02 ----
X <- cbind(1, as.matrix(casp[, -1]))
y <- casp$RMSD
QuanSub(X, y, r0 = 200, r = 1000, RR = 10, tau = 0.75, 
        optmethod = "L")





