2013-03-25 41 views
3

我有面板数据,并且在某些年份之前有很多变量缺少观测值。年份因变量而异。推断跨多列的缺失数据点的有效方法是什么?我正在考虑像线性趋势外推一样简单的事情,但我希望能找到一种将预测应用于多列的有效方法。下面是一个示例数据集,其缺失与我正在处理的类似。在这个例子中,我希望使用每列中观察到的数据点计算的线性趋势来填写“国家GDP”和“国民预期寿命”变量中的NA值。如何有效推断多个变量的缺失数据

###Simulate National GDP values 
set.seed(42) 
nat_gdp <- c(replicate(20L, { 
    foo <- rnorm(3, mean = 2000, sd = 300) + c(0,1000,2000) 
    c(NA,NA,foo)})) 
###Simulate national life expectancy values 

nat_life <- c(replicate(20L, { 
    foo <- rnorm(2, mean = 55, sd = 7.8) + c(0,1.5) 
    c(NA,NA,NA,foo)})) 




###Construct the data.table  
data.sim <- data.table( GovernorateID = c(rep(seq.int(11L,15L,by=1L), each = 20)), 
         DistrictID =rep(seq.int(1100,1500,by=100),each=20) + rep(seq_len(4), each = 5), 
         Year = seq.int(1990,1994,by=1L), 
         National_gdp = nat_gdp , 
         National_life_exp = nat_life ) 
+0

有你看着多重填补?参见“鼠标”包或其他相关软件包。此外,这不像统计问题那样是一个编程问题。 – ndoogan 2013-03-25 00:22:36

+0

永远不要在循环中增长矢量,预先分配或使用'replicate'。我编辑过,以显示更简单,更高效的方法 – mnel 2013-03-25 00:32:13

+0

@ndoogan多年来一直缺失,我错过了大多数变量的观察结果,因此多重插补不是一种选择。这就是为什么我在寻找简单的线性外推。 – 2013-03-25 00:34:38

回答

4

我假设你想分别对每个DistrictID做线性模型。

原始数据表:

head(data.sim) 
## GovernorateID DistrictID Year National_gdp National_life_exp 
## 1:   11  1101 1990   NA    NA 
## 2:   11  1101 1991   NA    NA 
## 3:   11  1101 1992  1988.746    NA 
## 4:   11  1101 1993  2527.619   54.70739 
## 5:   11  1101 1994  3854.210   44.21809 
## 6:   11  1102 1990   NA    NA 

dd <- copy(data.sim) # Make a copy for later. 

替换每个NA元素与线性模型的预测。一个链式操作中的两个步骤。

data.sim[, National_life_exp := ifelse(is.na(National_life_exp), 
             predict(lm(National_life_exp ~ Year, data=.SD), .SD), 
             National_life_exp) 
     , by=DistrictID 
     ][, National_gdp := ifelse(is.na(National_gdp), 
            predict(lm(National_gdp ~ Year, data=.SD), .SD), 
            National_gdp) 
      , by=DistrictID 
     ] 


head(data.sim) 
## GovernorateID DistrictID Year National_gdp National_life_exp 
## 1:   11  1101 1990 -8.004377   86.17531 
## 2:   11  1101 1991 924.727559   75.68601 
## 3:   11  1101 1992 1988.745871   65.19670 
## 4:   11  1101 1993 2527.618676   54.70739 
## 5:   11  1101 1994 3854.209743   44.21809 
## 6:   11  1102 1990 1008.886661   70.45643 

一个不错的(?)情节。请注意,在本例中,DistrictID的每个级别都有两个非NA点。

plot(data.sim$National_life_exp) 
points(dd$National_life_exp, col='red') # The copy from before. 

enter image description here

+0

这太棒了!谢谢!什么是“.SD”呢? – 2013-03-25 02:00:35

+0

这是用'by'选择'D'ata的'S'ubset。 – 2013-03-25 02:02:42