我正在尝试分析来自明尼苏达大学 IPUMS 数据集的1990年美国人口普查数据。R我正在使用该survey包，因为数据是加权的。仅获取家庭数据（并忽略人员变量以使事情简单），我正在尝试计算hhincome（家庭收入）的平均值。为此，我使用具有以下代码的函数创建了一个调查设计对象：svydesign()
> 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（使用用于同一数据集的不同部分的代码），我会得到不同的标准错误：
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了以下建议：
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
我似乎无法修复。这可能与这个问题有关。
总而言之：
- 为什么我在
Stataand中得不到相同的答案R？ - 哪一个是对的（或者我在这两种情况下都做错了什么）？
- 假设我得到了
rep()解决方案，那会复制Stata的结果吗？ - 正确的方法是什么？如果答案允许我使用
plyr包进行任意计算，而不是仅限于survey(等)中实现的函数svymean()，则表示敬意svyglm()
更新
因此，在我在这里和通过电子邮件从 IPUMS 收到的出色帮助之后，我正在使用以下代码来正确处理调查权重。我在这里描述以防其他人将来遇到此问题。
初始状态准备
由于 IPUMS 目前没有发布用于将数据导入 的脚本R，因此您需要从Stata、SAS或SPSS. 我暂时坚持Stata。首先从 IPUMS 运行导入脚本。然后在继续添加以下变量之前：
generate strata = statefip*100000 + puma
这将为 240001 表格中的每一个创建一个唯一的整数，PUMA前两位数字作为州 fip 代码（在马里兰州的情况下为 24），后四位PUMA是每个州唯一的 id。如果您要使用R，您可能还会发现运行它也很有帮助
generate statefip_num = statefip * 1
这将创建一个没有标签的附加变量，因为将.dta文件导入R应用标签并丢失底层整数。
统计数据和svyset
正如基思解释的那样，调查抽样是Stata通过调用来处理的svyset。
对于个人水平分析，我现在使用：
svyset serial [pweight=perwt], strata(strata)
这将权重设置为perwt，将分层设置为我们在上面创建的变量，并使用家庭serial人数来解释聚类。如果我们使用多年，我们可能想尝试
generate double yearserial = year*100000000 + serial
还要考虑纵向聚类。
对于家庭层面的分析（没有年份）：
svyset serial [pweight=hhwt], strata(strata)
应该是不言自明的（尽管我认为在这种情况下连续剧实际上是多余的）。替换serial为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。