library(foreign) #Reads Stata data library(MCMCpack) # An MCMC package for the ordered probit hrdata <- read.dta("~/data/hr/may2007/giant-dataset.dta") #Change the location to fit needs names(hrdata) # What variables do we have? ai.polity.base <- MCMCoprobit(aicombined ~ ailagv2+ailagv3+ailagv4+ailagv5+DEMOC+ai2democ+ai3democ+ai4democ+ai5democ+civwar01+intwar01+icvwar01+gdpgrowth+loggdppc+popchange+logpop, data=hrdata, mcmc=100000, burnin=25000) # This command line uses interaction terms that are already in the data. ai.polity.create <- MCMCoprobit(aicombined ~ ailagv2+ailagv3+ailagv4+ailagv5+DEMOC+as.factor(ailcomb):DEMOC+civwar01+intwar01+icvwar01+gdpgrowth+loggdppc+popchange+logpop, data=hrdata, mcmc=100000, burnin=25000) # This command line creates the interaction terms by treating lagged y as a factor and multiplying them by DEMOC. test.ai2DEMOC <- ai.polity.base[,"DEMOC"]+ai.polity.base[,"ai2democ"] # Create alpha(y(t-1)=2) + beta for Markov effects test.ai3DEMOC <- ai.polity.base[,"DEMOC"]+ai.polity.base[,"ai3democ"] # Create alpha(y(t-1)=3) + beta for Markov effects test.ai4DEMOC <- ai.polity.base[,"DEMOC"]+ai.polity.base[,"ai4democ"] # Create alpha(y(t-1)=4) + beta for Markov effects test.ai5DEMOC <- ai.polity.base[,"DEMOC"]+ai.polity.base[,"ai5democ"] # Create alpha(y(t-1)=5) + beta for Markov effects par(mfrow=c(2,2)) # Set up a graphics device to plot the four test densities plot(density(test.ai2DEMOC), main="AI Lag 2") plot(density(test.ai3DEMOC), main="AI Lag 3") plot(density(test.ai4DEMOC), main="AI Lag 4") plot(density(test.ai5DEMOC), main="AI Lag 5")