# Example 11.7, Wages and Productivity
# Data set: earns
# Function for result reporting
source("_report.R")
# Load the data and estimate the models in the background
load("earns.Rdata")
names(data)[names(data)=="lhrwage"]="lhrwage.t"
names(data)[names(data)=="loutphr"]="loutphr.t"
tsdata=ts(data,start=1947,frequency=1)
model1=lm(lhrwage.t~loutphr.t+t,data=tsdata)
model2=lm(ghrwage~goutphr,data=tsdata)
dig1=c(2,2,3,3)
dig2=c(4,3,3)
# Describe first model
cat("Model to estimate: lhrwage.t = beta0 + beta1 * loutphr.t + beta2 * t + u.t",
"\nwhere lhrwage is ", paste(desc[desc[,1]=="lhrwage",2]), " (hrwage: ", paste(desc[desc[,1]=="hrwage",2]), ")",
"\nand loutphr is ", paste(desc[desc[,1]=="loutphr",2]), " (outphr: ", paste(desc[desc[,1]=="outphr",2]), ")",
sep="")
# Report results
{
cat("The estimated regression line is")
reportreg(model1,dig1,suffix=".hat",adj=T)
}
# Interpretation
cat("The slope coefficient indicates that a 1% increase in outphr is predicted to raise hrwage by ",
printcoef(model1,2,dig1[2]), "%. Thus, the elasticity of hrwage respect to outphr suggested by the equation is quite large, and significantly different from unity",
"\nIn addition, we report the usual goodness-of-fit measures here, which do not correctly reflect the actual explaining power of the model. It would be better to report those based on the detrended dependent variable, as in Example 10.10",
sep="")
# Describe second model
m1=lm(lhrwage.t~t,data=tsdata) # These four lines obtain the
r1=ts(m1$residuals,start=1947,frequency=1) # autocorrelation of detrended
temp=as.data.frame(cbind(r1,r1_1=lag(r1,k=-1))) # lhrwage. Similar for the next
cor1=round(cor(temp$r1,temp$r1_1,use="complete.obs"),3) # four lines
m2=lm(loutphr.t~t,data=tsdata)
r2=ts(m2$residuals,start=1947,frequency=1)
temp=as.data.frame(cbind(r2,r2_1=lag(r2,k=-1)))
cor2=round(cor(temp$r2,temp$r2_1,use="complete.obs"),3)
cat("However, the first order autocorrelations of detrended lhrwage and loutphr are still as high as ",
cor1, " and ", cor2, ", respectively. To deal with the unit root problems, we estimate the equation again using first differences, ghrwage and goutphr (where g stands for \"growth\"), and without the time trend",
sep="")
# Report results
{
cat("The estimated regression line is")
reportreg(model2,dig2,adj=T)
}
# Interpretation
cat("This equation suggests that a 1% increase in outphr is predicted to raise hrwage by ",
printcoef(model2,2,dig2[2]), "%, and the estimate is not statistically different from unity. The adjusted R^2 indicates that the growth in outphr explains about ",
100*as.numeric(printrbar(model2,3)), "% of the growth in hrwage",
sep="")