我正在尝试分析明尼苏达大学 IPUMS 数据集中的数据1990 年美国人口普查 http://usa.ipums.org/usa/sampdesc.shtml#us1990a in R。我正在使用survey http://faculty.washington.edu/tlumley/survey/包，因为数据是weighted http://en.wikipedia.org/wiki/Sampling_%28statistics%29#Survey_weights。仅获取家庭数据（并忽略个人变量以使事情简单），我试图计算hhincome (家庭收入 http://internal.usa.ipums.org/usa-action/variables/HHINCOME）。为此，我创建了一个调查设计对象使用svydesign() http://faculty.washington.edu/tlumley/survey/example-design.html函数具有以下代码：
> require(foreign)
> ipums.household <- read.dta("/path/to/stata_export.dta")
> ipums.household[ipums.household$hhincome==9999999, "hhincome"] <- NA # Fix missing
> ipums.hh.design <- svydesign(id=~1, weights=~hhwt, data=ipums.household)
> svymean(ipums.household$hhincome, ipums.hh.design, na.rm=TRUE)
mean SE
[1,] 37029 17.365
到目前为止，一切都很好。但是，如果我尝试在中进行相同的计算，我会得到不同的标准错误Stata (using 用于同一数据集不同部分的代码 http://www.stanford.edu/~mrosenfe/soc_meth_proj3/Intro%20to%20STATA%20for%20Soc%20180.htm):
use "C:\I\Hate\Backslashes\stata_export.dta"
replace hhincome = . if hhincome == 9999999
(933734 real changes made, 933734 to missing)
mean hhincome [fweight = hhwt] # The code from the link above.
Mean estimation Number of obs = 91746420
--------------------------------------------------------------
| Mean Std. Err. [95% Conf. Interval]
-------------+------------------------------------------------
hhincome | 37028.99 3.542749 37022.05 37035.94
--------------------------------------------------------------
而且，在寻找另一种给这只猫剥皮的方法时，作者survey, has 这个建议 http://pj.freefaculty.org/R/Rtips.html#1.12对于频率加权：
expanded.data<-as.data.frame(lapply(compressed.data,
function(x) rep(x,compressed.data$weights)))
但是，我似乎无法让这段代码工作：
> hh.dataframe <- data.frame(ipums.household$hhincome, ipums.household$hhwt)
> expanded.hh.dataframe <- as.data.frame(lapply(hh.dataframe, function(x) rep(x, hh.dataframe$hhwt)))
Error in rep(x, hh.dataframe$hhwt) : invalid 'times' argument
我似乎无法解决这个问题。这可能与这个问题 http://r.789695.n4.nabble.com/error-with-source-invalid-times-value-td3234425.html.
总而言之：
- 为什么我没有得到相同的答案
Stata and R?
- 哪一个是正确的（或者我在这两种情况下都做错了什么）？
- 假设我得到了
rep()解决方案有效，会复制吗Stata的结果？
- 什么是right怎么办？如果答案允许我使用，那就太好了plyr http://had.co.nz/plyr/用于进行任意计算的包，而不是仅限于中实现的功能
survey (svymean(), svyglm() etc.)
Update
因此，在我通过电子邮件从 IPUMS 获得了出色的帮助后，我使用以下代码来正确处理调查权重。我在这里描述以防其他人将来遇到这个问题。
初始 Stata 准备
由于 IPUMS 目前不发布用于将数据导入到的脚本R，你需要从Stata, SAS, or SPSS。我会坚持Stata目前。首先从 IPUMS 运行导入脚本。然后在继续之前添加以下变量：
generate strata = statefip*100000 + puma
这为每个创建一个唯一的整数PUMA格式为 240001，前两位数字为州 FIP 代码（马里兰州为 24），后四位数字为PUMAid 在每个州都是唯一的。如果你要使用R您可能还会发现运行它也很有帮助
generate statefip_num = statefip * 1
这将创建一个没有标签的附加变量，因为导入.dta文件到R应用标签并丢失底层整数。
斯塔塔和svyset
正如 Keith 所解释的，调查抽样是由Stata通过调用svyset.
对于个人层面的分析，我现在使用：
svyset serial [pweight=perwt], strata(strata)
这将权重设置为perwt，对我们上面创建的变量进行分层，并使用家庭serial用于说明聚类的数字。如果我们使用多年，我们可能想尝试
generate double yearserial = year*100000000 + serial
也考虑纵向聚类。
对于家庭层面的分析（不包括年份）：
svyset serial [pweight=hhwt], strata(strata)
应该是不言自明的（尽管我认为在这种情况下序列实际上是多余的）。更换serial with yearserial将考虑时间序列。
这样做在R
假设您要导入.dta附加文件strata上面解释的变量并在单个字母处进行分析：
require(foreign)
ipums <- read.dta('/path/to/data.dta')
require(survey)
ipums.design <- svydesign(id=~serial, strata=~strata, data=ipums, weights=perwt)
或者在家庭层面：
ipums.hh.design <- svydesign(id=~serial, strata=~strata, data=ipums, weights=hhwt)
希望有人觉得这有帮助，并非常感谢 IPUMS 的 Dwin、Keith 和 Brandon。
