R - 回归中的单位特定时间趋势

Question

在回归中，我试图对单位特定的时间趋势进行建模，但我一直遇到困难。在R当我估计用单位和年份固定效应模型就像lm(y~x+factor(unit)+factor(time))我得到完全正常的结果。但是，当我尝试这样做时，lm(y~x+factor(unit)*factor(year))我遇到了NA's所产生的麻烦。

使用一些模拟数据来说明：

# Unit of analysis are countries
country<-c(rep("Isthmus",10),rep("Nambutu",10),rep("San Monique",10))
ccode<-c(rep(1,10),rep(2,10),rep(3,10))
year <- c(rep(2000:2009,3)) # Time
x1<-rnorm(30)*ccode 
x2<-runif(30)
y<-0.5*x1-0.3*x2+rnorm(30) # Outcome variable
df=data.frame(country,ccode,year,y,x1,x2)

分别使用单位和时间、国家和年份的固定效应来估计模型：

m0<-lm(y~x1+x2+factor(ccode)+factor(year),df);summary(m0)

# Part of the regression output:

Coefficients:
                  Estimate Std. Error t value  Pr(>|t|)    
(Intercept)      -0.92780    0.68231  -1.360    0.1928    
x1                0.59290    0.10058   5.895 0.0000226 ***
x2               -0.36457    0.96036  -0.380    0.7092    
factor(ccode)2    0.95383    0.48675   1.960    0.0677 .  
factor(ccode)3    0.46050    0.46475   0.991    0.3365    
factor(year)2001  0.15222    0.87295   0.174    0.8638 

这里没有问题。现在我使用特定于单位的时间趋势来估计模型：

m1<-lm(y~x1+x2+factor(year)*factor(ccode),df);summary(m1)

# Part of the regression output:

                                  Estimate Std. Error t value Pr(>|t|)
(Intercept)                       1.3408         NA      NA       NA
x1                                3.3104         NA      NA       NA
x2                                0.5239         NA      NA       NA
factor(ccode)2                   -2.0544         NA      NA       NA
factor(ccode)3                  -12.2971         NA      NA       NA
factor(year)2001:factor(ccode)1  -3.4409         NA      NA       NA
factor(year)2002:factor(ccode)1  -0.6348         NA      NA       NA

在这种特殊情况下，这NA's似乎是模型中变量过多的结果，因为没有剩余的自由度。使用更大的数据集时也会出现同样的问题。我不完全确定这里出了什么问题。我认为这与我factor用来对单位特定时间趋势进行建模的方式有关，但到目前为止我还没有解决这个问题。

有没有人知道如何正确地做到这一点？欢迎任何建议。

Answer 1

您正在尝试估计比数据更多的参数，即n < p。在您的示例数据集中，您有

R> nrow(df)
[1] 30

数据点，并试图估计 30 个参数。正如 Ben 指出的那样，您估计每年都有不同的参数。如果你想假设一个线性趋势，那么只要

lm(y ~ x1 + factor(ccode)*time, data=df)

或包含二次趋势

lm(y ~ x1 + factor(ccode)*I(time^2), data=df)