# \donttest{
#Example 1
##Load the data
data(ais)
attach(ais)
## Fitting a median regression with Normal errors (by default)
modelF = lqr(BMI~LBM,data = ais,subset = "(Sex==1)")
modelM = lqr(BMI~LBM,data = ais,subset = "(Sex==0)")
plot(LBM,BMI,col=Sex*2+1,
xlab="Lean Body Mass",
ylab="Body4 Mass Index",
main="Quantile Regression")
abline(a = modelF$beta[1],b = modelF$beta[2],lwd=2,col=3)
abline(a = modelM$beta[1],b = modelM$beta[2],lwd=2,col=1)
legend(x = "topleft",legend = c("Male","Female"),lwd = 2,col = c(1,3))
#COMPARING SOME MODELS for median regression
modelN = lqr(BMI~LBM,dist = "normal")
modelT = lqr(BMI~LBM,dist = "t")
modelL = lqr(BMI~LBM,dist = "laplace")
#Comparing AIC criteria
modelN$AIC;modelT$AIC;modelL$AIC
#This could be automatically done using best.lqr()
best.model = best.lqr(BMI~LBM,data = ais,
p = 0.75, #third quartile
criterion = "AIC")
#Let's use a grid of quantiles (no output)
modelfull = lqr(BMI~LBM,data = ais,
p = seq(from = 0.10,to = 0.90,by = 0.05),
dist = "normal",silent = TRUE)
#Plotting quantiles 0.10,0.25,0.50,0.75 and 0.90
if(TRUE){
plot(LBM,BMI,xlab = "Lean Body Mass"
,ylab = "Body Mass Index", main = "Quantile Regression",pch=16)
colvec = c(2,2,3,3,4)
imodel = c(1,17,4,14,9)
for(i in 1:5){
abline(a = modelfull[[imodel[i]]]$beta[1],
b = modelfull[[imodel[i]]]$beta[2],
lwd=2,col=colvec[i])
}
legend(x = "topleft",
legend = rev(c("0.10","0.25","0.50","0.75","0.90")),
lwd = 2,col = c(2,3,4,3,2))
}
#Example 2
##Load the data
data(crabs,package = "MASS")
attach(crabs)
## Fitting a median regression with Normal errors (by default) #Note the double quotes
crabsF = lqr(BD~FL,data = crabs,subset = "(sex=='F')")
crabsM = lqr(BD~FL,data = crabs,subset = "(sex=='M')")
if(TRUE){
plot(FL,BD,col=as.numeric(sex)+1,
xlab="Frontal lobe size",ylab="Body depth",main="Quantile Regression")
abline(a = crabsF$beta[1],b = crabsF$beta[2],lwd=2,col=2)
abline(a = crabsM$beta[1],b = crabsM$beta[2],lwd=2,col=3)
legend(x = "topleft",legend = c("Male","Female"),
lwd = 2,col = c(3,2))
}
#Median regression for different distributions
modelN = lqr(BD~FL,dist = "normal")
modelT = lqr(BD~FL,dist = "t")
modelL = lqr(BD~FL,dist = "laplace")
modelS = lqr(BD~FL,dist = "slash")
modelC = lqr(BD~FL,dist = "cont" )
#Comparing AIC criterias
modelN$AIC;modelT$AIC;modelL$AIC;modelS$AIC;modelC$AIC
# best model based on BIC
best.lqr(BD~FL,criterion = "BIC")
#Let's use a grid of quantiles for the Student's t distribution
modelfull = lqr(BD~FL,data = crabs,
p = seq(from = 0.10,to = 0.90,by = 0.05),
dist = "t") # silent = FALSE
#Plotting quantiles 0.10,0.25,0.50,0.75 and 0.90
if(TRUE){
plot(FL,BD,xlab = "Frontal lobe size"
,ylab = "Body depth", main = "Quantile Regression",pch=16)
colvec = c(2,2,3,3,4)
imodel = c(1,17,4,14,9)
for(i in 1:5){
abline(a = modelfull[[imodel[i]]]$beta[1],
b = modelfull[[imodel[i]]]$beta[2],
lwd=2,col=colvec[i])
}
legend(x = "topleft",
legend = rev(c("0.10","0.25","0.50","0.75","0.90")),
lwd = 2,col = c(2,3,4,3,2))
}
# }
Run the code above in your browser using DataLab