现在大家都知道我强烈推荐在DID中进行“平行趋势敏感性分析”。这是因为Roth（2022，AER：Insights）指出，传统的处理前趋势检验效力在50%-80%之间。这个效力是什么？这就是今天想给大家介绍的，也是我的第三篇方法论模拟分析文章《因果推断的效力分析的zui新趋势与实践指导》（前两篇分别是《听Rubin和Imbens的话：处理配置机制讨论的实践指导》和《如何控制协变量？DID中协变量调整的实践指导》）。
除了Roth的文章外，zui近几年，效力分析也在Top期刊的应用研究中越来越常见，例如，Black et al.（2022，JPubE）、Dench et al.（2024，JPubE）。关于效力分析更详细的技术讲解，请参见Sylvain Chabé-Ferret（2024）的计量课本《Statistical Tools for Causal Inference》第七章。
下面的内容主要来自于Nick Huntington-Klein（2020）：《Simulation for Power Analysis》。
我们有了idea后，就要去做许多统计分析。在学习计量之初，老师就告诉我们，要尽可能的增加样本量。效力分析（Power Analysis）就是一种我们至少需要多少样本量得到给定处理效应，或者给定样本量可以得到zui低处理效应是多少。
无论你进行哪种研究，效力分析都是一个好主意。然而，它在两种情况下特别有用：
如果你正在寻找一个你认为可能不是那么核心或重要的影响，即它是一个小的影响，或是一个系统的一部分，在这个系统中，有很多其他的事情正在发生（粗略地，你可以把它想象成“一个小的 2 ”），进行效力分析可能是一个好主意——了解小效应所需的样本量通常比你预期的要大得多，zui好现在就了解这一点，而不是在你已经完成所有工作之后再做
如果你正在进行随机实验，你实际上可以控制样本量 - 你可以选择收集多少数据以及如何随机化处理。在进行实验之前，效力分析至关重要，这样你就不会在实验结束时才意识到“哦，该死，我应该再多几百人做一次......现在太晚了！”
效力分析的作用
让X表示处理，Y表示结果，假设我们正在研究X和Y关系，效力分析设计以下五件事：
① 处理效应的大小
② 处理的变动
③ 其它因素引起Y的变动
④ 统计精度
⑤ 样本量
效力分析就是保持其它四项不变，第五项是多少。例如，如果我们认为效应可能是 A，并且X的变动为B，并且Y中有 C 的变动与X无关，并且你希望至少有 D% 的机会发现一个效应（如果确实存在），那么你需要至少 E 的样本量。这告诉我们获得足够统计效力所需的zui小样本量是E。
或者如果你有一个样本量为 A，并且X的变动为 B ，并且Y有 C 的变动与X无关，并且你希望至少有 D% 的机会发现一个效应（如果它确实存在），那么这个效应必须至少与 E 一样大。这告诉我们可检测到的效应的下界，即在给定样本大小的情况下，你希望有机会合理测量的zui小效应。
那么，“统计精度”是什么呢？通常，你有一个目标统计效力水平（因此称为“效力分析”）。统计效力是真实率。也就是说，如果确实存在影响，并且抽样变动意味着你有 80% 的机会拒绝给定样本中无效的零假设，那么你就有 80% 的统计效力。统计效力还取决于你正在运行的测试类型 - 如果你在 95% 的置信水平下进行假设检验，则与在 99% 的置信水平下进行假设检验相比，你更有可能拒绝零假设（因此统计效力更高）。
效力分析并不一定要考虑统计效力。我们可以从以上五项的任何一个角度来进行效力分析，例如标准误。
这五项数值的来源
为了进行效力分析，我们需要知道五个数值中的四个，这样效力分析才能告诉我们第五个。
另一个问题是效力标准是多大？统计效力越高，我们错过真实处理效应的可能性越小，但这也意味着我们需要更大的样本。一般的经验来看，80%的效力是一个标准，但现在越来越多的研究者使用90%的效力。
效力分析的模拟
第yi步：构造模拟数据集
我们模拟1000个样本，处理变量X是0,1上的均匀分布，处理效应是0.2，标准误是0-3的正态分布。stata代码如下：
* Set the number of observations we want in our simulated dataclear* Create our data* Note 0 and 1 are the default start and end of runiform; I want that, so I don't put anything inside runiform()g X = runiform()* Now create Y based on X* Don't forget to add additional noise!* The *true effect* of X on Y is .2g Y = .2*X + rnormal(0, 3)
第二步：执行分析
假设我们计划获得回归中的稳健标准误：
* reg is short for regressreg Y X, robust* And if we're just interested in significance, say at the 95% level,* we can compare the p-value to .05 and store the result as a local variable (i.e. just store the single number, not as a regular Stata variable)local sig = 2*ttail(e(df_r),abs(_b[X]/_se[X])) <= .05* di is short for displaydi `sig'
第三步：重复上述过程
forvalues i = 1(1)500 { * Set the number of observations we want in our simulated data clear set obs 1000 * Create our data g X = runiform() g Y = .2*X + rnormal(0, 3) * Run analysis quietly qui reg Y X, robust * Pull out results local sig = 2*ttail(e(df_r),abs(_b[X]/_se[X])) <= .05 di `sig' local coef = _b[X] di `coef'}这个重复过程将会产生500次数据，然后他将获取X的系数以及95%置信水平上是否显著。
第四步：储存结果
我们首先创建一个数据集来存储结果。然后，随着我们反复进行分析，我们将该preserve结果数据集清除，进行数据生成，并将结果存储在本地变量中。然后我们restore返回结果数据集并将这些本地变量放入实际变量中。zui后，我们将拥有一个充满结果的数据集！
clear* We want a number of observations equal to the number of times we will simulate the dataset obs 500* Blank variablesg coef_results = .g sig_results = .* Let's wrap the whole thing in quietly{} - we don't need the output!quietly { forvalues i = 1(1)500 { * Preserve our results data set so we can instantly come back to it preserve * Set the number of observations we want in our simulated data * Now we're no longer in the results data set; we're in the simulated-data data set clear set obs 1000 * Create our data g X = runiform() g Y = .2*X + rnormal(0, 3) * Run analysis quietly qui reg Y X, robust * Pull out results local sig = 2*ttail(e(df_r),abs(_b[X]/_se[X])) <= .05 local coef = _b[X] * Now restore to come back to results restore * And store our results in the results variables, using "in" to tell it which row to store the data in replace coef_results = `coef' in `i' replace sig_results = `sig' in `i' }}* summ is short for summarizesumm sig_results结果如下：
结果显示，6%的统计效力（均值）。
第五步：更多的结果
我们还可以考察系数的分布：
tw kdensity coef_results, xti("Coefficient on X") yti("Density") lc(red)X系数的密度分布
label def sig 0 "Not Significant" 1 "Significant"label values sig_results siggraph bar, over(sig_results) yti("Percentage")第七步：高级效力分析
效力分析还可以获得更多的信息：（1）我们想要分析不同效应下的效力是多少。首先固定1000个样本，尝试分析不同效应的大小的效力水平：
foreach effect in .4 .8 1.2 1.6 2 { foreach sample_size in 1000 { quietly { clear * We want a number of observations equal to the number of times we will simulate the data set obs 500 * Blank variables g sig_results = . * Let's wrap the whole thing in quietly{} - we don't need the output! forvalues i = 1(1)500 { preserve clear * NOTICE I'm putting the "sample_size" variable from above here set obs `sample_size' * Create our data g X = runiform() * and the "effect" variable here g Y = `effect'*X + rnormal(0, 3) * Run analysis quietly reg Y X, robust * Pull out results local sig = 2*ttail(e(df_r),abs(_b[X]/_se[X])) <= .05 local coef = _b[X] * Now restore to come back to results restore * And store our results in the results variables, using "in" to tell it which row to store the data in replace sig_results = `sig' in `i' } summ sig_results * since we're inside a quietly{} we need a noisily to see this noisily di "With a sample size of `sample_size' and an effect of `effect', the mean of sig_results is " r(mean) } }}因此，看起来我们需要一个介于 0.8 和 1.2 之间的效应，才能有 80% 的机会找到显著的结果。如果我们认为效应实际上不太可能那么大，那么我们需要想其他办法---找到更大的样本，使用更精确的估计方法，等等！否则我们可能应该放弃。
