Multiple Regression

#############################
#BLOCK #1
library('BayesFactor')
#load data
#dat=read.table('bailey.dat',head=T)
dat=read.table(url("http://pcl.missouri.edu/jeff/sites/pcl.missouri.edu.jeff/files/bailey.dat"),head=T)
#############################

#############################
#BLOCK #2, Conventional Analysis

cc=dat$CC
climate=dat$diftemp
pop=dat$dpop_30
para=dat$tpar
o2=dat$Isosd

#make a new DataFrame, needed for BF
myDat=data.frame(cc,climate,pop,para,o2)

summary(lm(cc~.,data=myDat))

round(cor(myDat),2)

#######################################
#BLOCK 3

bfTop=regressionBF(cc~climate+para+o2+pop,data=myDat,whichModels="top")

bfTop

plot(bfTop)

#######################################
#BLOCK 4 Bayes Factor, All Models, Go Slow

bf=regressionBF(cc~climate+para+o2+pop,data=myDat)

length(bf)
head(bf)
sort(bf)
tail(bf)
plot(bf)

bf[11]
bf[11:15]

#renormalizng BF
bf

bf.normalized=bf/bf[15]

plot(bf.normalized)

bf[10]/bf[14] #better than the addition of para
bf[10]/bf[13] #better than the addition of climate
bf[10]/bf[3] #better than the omission of pop
bf[10]/bf[4] #better than the omission of 02