手动计算 z 检验的功效二元变量

问题描述

我正在尝试计算最近发表的一项研究的功效,其中进行了 z 检验以证明免疫抑制药物在感染 sars-cov-2 后显着增加了初步病死率。 但是,来自 statsmodels.stats.power 的 zt_ind_solve_power 给出了不同的结果,我不明白为什么! 这就是我所做的!我根据组大小和标准偏差计算了总体中的预期标准误差,然后从临界值(此处:临界值基于抽取的样本)计算了替代曲线下的面积(此处:近似分布,基于抽取的样本)论文中报道了p):

import numpy as np
import scipy.stats as stats
import scipy.integrate as integrate
from statsmodels.stats.power import zt_ind_solve_power
n1 = 29
n2 = 37

mean1 = 0.28 # mortaility rate in control
mean2 = 0.7 # mortaility rate in treatment

x1 = n1*mean1
x2 = n2*mean2

p = (x1+x2)/(n1+n2)

sd1 = n1*P*(1-p)
sd2 = n2*P*(1-p)

sdp = np.sqrt(((n1-1)*sd1**2+(n2-1)*sd2**2)/n1+n2-2)
d = (mean2-mean1)/sdp

mean_null = 0 # H0: mortality rate equal in both groups
mean_sample = mean2-mean1

se = np.sqrt(P*(1-p)*((1/n1)+(1/n2)))
x = np.linspace(-0.5,1,num=10000)
null_dist = stats.norm.pdf(x=x,loc=mean_null,scale=se)
sample_dist = stats.norm.pdf(x=x,loc=mean_sample,scale=se)
p_val = 0.0013
critical_one = stats.norm.sf(x=p_val,scale=se)# assuming one-tailed test!

for i in range(x.size):
    if x[i] > critical_one and x[i-1] <= critical_one:
        print(x[i]-critical_one)
        index_one = i
    if x[i] > critical_two and x[i-1] <= critical_two:
        index_two = i
        print(x[i]-critical_two)
power_one = integrate.simps(y=sample_dist[index_one:x.size+1],x=x[index_one:x.size+1])

作为我的学士论文的一部分,我正在尝试手动求解权力(我正在分析各种研究人员的自由度)。我知道,否则,没有必要这样做,真的…… 让我想知道的是,solve_power 函数在 p=0.0013 和 power=6.6% 时返回 power=1.9%,这似乎太低了,在这里不真实。 在 p=1.3% 时,我的功率=27%。

zt_ind_solve_power(effect_size = d,nobs1 = n1,alpha=p_val,ratio=n2/n1,alternative='larger')

我无法理解我做错了什么!我不希望 zt_ind_solve_power 犯错,所以...你能告诉我我错过了什么吗?

非常感谢!!!

解决方法

暂无找到可以解决该程序问题的有效方法,小编努力寻找整理中!

如果你已经找到好的解决方法,欢迎将解决方案带上本链接一起发送给小编。

小编邮箱:dio#foxmail.com (将#修改为@)