counter<-function(formel,data){ ###################################################################################### library(quantreg) sample<-data # SET TAUS taus<-seq(0.05,0.95,0.05) # ESTIMATE QR qreg_0<-rq(formula=formel,data=subset(sample,group==0),tau=taus) qreg_1<-rq(formula=formel,data=subset(sample,group==1),tau=taus) #summary(qreg_0) #summary(qreg_1) # STORE PREDICTED VALUES FOR EACH TAU IN EACH GROUP pred0 <- data.frame(predict(qreg_0,sample)) names(pred0) <- gsub("tau\\.\\.0\\.", "q_0_", names(pred0)) names(pred0) <- gsub("q\\_0\\_0", "q\\_0\\_", names(pred0)) pred1 <- data.frame(predict(qreg_1,sample)) names(pred1) <- gsub("tau\\.\\.0\\.", "q_1_", names(pred1)) names(pred1) <- gsub("q\\_1\\_0", "q\\_1\\_", names(pred1)) # COMPBINE ESTIMATED VALUES OF GROUP 0 AND 1 prediction<-cbind(group=sample$group,pred0,pred1) #summary(prediction) prediction[is.na(prediction$group),]<-NA # PREPARE TAUS ctaus<-as.character(taus*100) # CALCULATE RELATIVE CHANGES FOR EACH TAU for(i in ctaus){ prediction[[paste('delta',i,sep="_")]]<-prediction[[paste('q_0',i,sep="_")]]/prediction[[paste('q_1',i,sep="_")]]} ###################################################################################### prediction$mieteqm<-(sample$mieteqm) prediction[is.na(prediction$group),]<-NA prediction <- within(prediction, tau <- 5*(as.integer(cut(mieteqm, quantile(mieteqm, probs=0:19/19,na.rm=T), include.lowest=TRUE)))) for(i in ctaus){ prediction[[paste('counter',i,sep="_")]]<-prediction[[paste('delta',i,sep="_")]]*prediction$mieteqm} ###################################################################################### prediction$counter<-NA for (i in 1:nrow(prediction)){ for(ct in ctaus){ if (prediction$tau[i]==ct & !is.na(prediction$group[i])) prediction$counter[i]<-prediction[[paste('counter',ct,sep="_")]][i] } } counterresults<<-cbind(data,counter=ifelse(prediction$group==1,prediction$counter,NA)) ###################################################################################### ###################################################################################### ###################################################################################### }