#reads in the data and splits it into various data frames alldata<-read.table("c:\\Myfiles\\Personal\\Papers\\RegSubset\\alldata.txt",header=T) base<-alldata[1:12,] train<-alldata[1:7,] predc<-alldata[8:12,] proj<-alldata[13:19,] #gets default graphic parameters defpar<-par() #fits the 5 variable model to the base period m1<-lm(EnergyC~X1+X2+X3+X4+X5, data=base) summary(m1) #makes plot of actual and predicted values for 5 variable model in base period ivv<-rep(c(5000,25000),c(6,6)) #modifies plot margins par(mai=c(2.5,.75,.5,.25)) plot(1:12,ivv,type="n",xlab="Month",ylab="Energy Consumption") axis(1,1:12) points(1:12,base$EnergyC,type="p",pch=1) lines(1:12,fitted(m1), lty=2) #fits the full model to the data in the training set m1<-lm(EnergyC~X1+X2+X3+X4+X5,train) summary(m1) #gets predicted values and residuals from full model in validation set p1p<-predict(m1,newdata=predc) r1p<-predc$EnergyC-p1p res1<-c(r1p,residuals(m1)) #gets predicted values from model 1 in period after improvements pr1p<-predict(m1,newdata=proj) #fits 3 variable model to data in training set m2<-lm(EnergyC~X1+X2+X5,train) summary(m2) #gets predicted values and residuals from 3 variable model in validation set p2p<-predict(m2,newdata=predc) r2p<-predc$EnergyC-p2p res2<-c(r2p,residuals(m2)) #gets predicted values from model 2 in project period after improvements pr2p<-predict(m2,newdata=proj) #restores plot margins to default par(defpar) #sets plot to have two plots in one column for Figure 3 par(mfrow=c(2,1)) #creates plot for 3 variable model #creates artificial variable to set axis ivv<-rep(c(8000,28000),c(12,7)) #makes plot with no points in it plot(1:19,ivv,type="n", xlab=" Training Validation After Improvements", ylab="Energy Consumed", main="Predictions from 3-variable model") axis(1,1:19) lines(c(7.5,7.5),c(8000,28000),lty=2) lines(c(12.5,12.5),c(8000,28000),lty=2) #adds actual points during training period points(1:7,train$EnergyC, type="p", pch=15) #adds actual points in prediction period points(8:12,predc$EnergyC, type="p", pch=17) #adds actual points during the project period points(13:19,proj$EnergyC, type="p", pch=18) #adds predicted values in training period lines(1:7,fitted(m2), lty=1) #adds predicted values in predc period lines(8:12,p2p, lty=1) #adds predicted values in project period lines(13:19,pr2p, lty=1) #creates plot for full model ivv<-rep(c(8000,28000),c(12,7)) #makes plot with no points in it plot(1:19,ivv,type="n", xlab=" Training Validation After Improvements", ylab="Energy Consumed",main="Predictions from 5-variable model") axis(1,1:19) text(9.5,-5000,"Month") lines(c(7.5,7.5),c(8000,28000),lty=2) lines(c(12.5,12.5),c(8000,28000),lty=2) #adds actual points during base period points(1:7,train$EnergyC, type="p", pch=15) #adds actual points in prediction period points(8:12,predc$EnergyC, type="p", pch=17) #adds actual points during the project period points(13:19,proj$EnergyC, type="p", pch=18) #adds predicted values in Base period lines(1:7,fitted(m1), lty=1) #adds predicted values in prediction period lines(8:12,p1p, lty=1) #adds predicted values in project period lines(13:19,pr1p, lty=1) #sets plot back to one plot for Figure 2 par(defpar) #par(mfrow=c(1,1)) #creates factor variable for grouping data ind<-factor(rep(c("5 Variable Model","3 Variable Model"),c(5,5))) #combines residuals from full model and 3 variable subset res<-c(r1p,r2p) #makes box plots of actual minus predicted in validation set plot(ind,res,ylab="Difference in Actual and Predicted", main="Precision of Prediction") #fits the full model to the data in the base period m1<-lm(EnergyC~X1+X2+X3+X4+X5,base) summary(m1) #gets predicted values from model 1 in period after improvements pr1p<-predict(m1,newdata=proj) #fits 3 variable model to data in base period m2<-lm(EnergyC~X1+X2+X5,base) summary(m2) #gets predicted values from model 2 in period after improvements pr2p<-predict(m2,newdata=proj) #sets plot to have two plots in one column for Figure 4 par(mfrow=c(2,1)) #creates plot for 3 variable model #creates artificial variable to set axis ivv<-rep(c(8000,18000),c(12,7)) #makes plot with no points in it plot(1:19,ivv,type="n", xlab=" Base Period After Improvements", ylab="Energy Consumed", main="Predictions from 3-variable model") axis(1,1:19) lines(c(12.5,12.5),c(8000,18000),lty=2) #adds actual points during training period points(1:7,train$EnergyC, type="p", pch=15) #adds actual points in prediction period points(8:12,predc$EnergyC, type="p", pch=17) #adds actual points during the project period points(13:19,proj$EnergyC, type="p", pch=18) #adds predicted values in base period lines(1:12,fitted(m2), lty=1) #adds predicted values in project period lines(13:19,pr2p, lty=1) #creates plot for full model ivv<-rep(c(8000,18000),c(12,7)) #makes plot with no points in it plot(1:19,ivv,type="n", xlab=" Base Period After Improvements", ylab="Energy Consumed",main="Predictions from 5-variable model") axis(1,1:19) lines(c(12.5,12.5),c(8000,18000),lty=2) #adds actual points during base period points(1:7,train$EnergyC, type="p", pch=15) #adds actual points in prediction period points(8:12,predc$EnergyC, type="p", pch=17) #adds actual points during the project period points(13:19,proj$EnergyC, type="p", pch=18) #adds predicted values in Base period lines(1:12,fitted(m1), lty=1) #adds predicted values in project period lines(13:19,pr1p, lty=1) #calculates predicted energy savings in period after improvements for each model rr1p<-pr1p-proj$EnergyC rr2p<-pr2p-proj$EnergyC cbind(13:19,proj$EnergyC,pr1p,pr2p,rr1p,rr2p) # calculates predicted cumulative energy savings for each model sum(rr1p)