# Example 4.5, Housing Prices and Air Pollution
# Data set: hprice2
# Function for result reporting
source("_report.R")
# Load the data and estimate the model in the background
load("hprice2.Rdata")
data$ldist=log(data$dist) # There is a variable missing in the
desc=rbind(desc,data.frame(variable="ldist",label="log(dist)")) # data set. These lines create it.
model=lm(lprice~lnox+ldist+rooms+stratio, data=data)
dig=c(2,3,3,3,3,3)
# Describe the model
cat("Model to estimate: lprice = beta0 + beta1 * lnox + beta2 * ldist + beta3 * rooms + beta4 * stratio + u",
"\nwhere lprice is ", paste(desc[desc[,1]=="lprice",2]), " (price: ", paste(desc[desc[,1]=="price",2]), ")",
"\nlnox is ", paste(desc[desc[,1]=="lnox",2]), " (nox: ", paste(desc[desc[,1]=="nox",2]), ")",
"\nldist is ", paste(desc[desc[,1]=="ldist",2]), " (dist: ", paste(desc[desc[,1]=="dist",2]), ")",
"\nrooms is ", paste(desc[desc[,1]=="rooms",2]),
"\nand stratio is ", paste(desc[desc[,1]=="stratio",2]),
sep="")
# Report results
{
cat("The estimated regression line is")
reportreg(model,dig)
}
# Interpretation
cat("The fact that beta1hat < 0 implies that price is predict to fall as nox rises. Since beta1 is the elasticity of price with respect to nox, we do not test H0: beta1 = 0, but test H0: beta1 = -1 against H1: beta1 ≠ -1",
"\nTo test H0: beta1 = -1, we compute the t statistic as (",
printcoef(model,2,dig[2]), " + 1)/", printse(model,2,dig[2]), " = ",
round((as.numeric(printcoef(model,2,dig[2]))+1)/as.numeric(printse(model,2,dig[2])),2),
", which is a very small number. Therefore, we fail to reject H0 even at very large significance levels, i.e. there is little evidence that the elasticity is different from unity",
sep="")