發表文章

目前顯示的是 1月, 2024的文章

Collider simulation

 # Setting seed for reproducibility set.seed(123) # Number of observations n = 1000 # Generating independent variables X1 and X2 X1 = rnorm(n) X2 = rnorm(n) # Positive Correlation Scenario # Collider positively influenced by both X1 and X2 Collider_Positive = X1 + X2 + rnorm(n) # Negative Correlation Scenario # Collider positively influenced by X1 and negatively by X2 Collider_Negative = X1 - X2 + rnorm(n) # Checking correlations # Before conditioning cat("Correlation between X1 and X2 before conditioning: ", cor(X1, X2), "\n") # After conditioning - Positive Scenario cat("Positive Scenario - After conditioning on the collider: ",     cor(X1[Collider_Positive > median(Collider_Positive)],          X2[Collider_Positive > median(Collider_Positive)]), "\n") # After conditioning - Negative Scenario cat("Negative Scenario - After conditioning on the collider: ",     cor(X1[Collider_Negative > median(Collider_Negative)],    ...

用 Python 模擬贏家的詛咒

 # 用 Python 模擬贏家的詛咒 # 贏家的詛咒:低訊號雜訊比、低效果量、低統計檢定力、小樣本的研究容易高估效果量 # 訊號雜訊比: 實際效果量/標準誤 import numpy as np import matplotlib.pyplot as plt import seaborn as sns from scipy.stats import ttest_ind # 模擬設定 np.random.seed(0) num_experiments = 1000  # 實驗次數 sample_size = 30       # 樣本大小 true_effect_size = 0.3  # 實際效果量(小) signal_to_noise_ratio = 0.5  # 訊號雜訊比(低) # 模擬實驗的函數 def simulate_experiment(): # 為控制組和治療組生成數據     control = np.random.normal(0, 1/signal_to_noise_ratio, sample_size)     treatment = np.random.normal(true_effect_size, 1/signal_to_noise_ratio, sample_size) # 進行t檢定 t_stat, p_val = ttest_ind(treatment, control) return p_val, np.mean(treatment) - np.mean(control) # 運行模擬 p_values = [] effect_sizes = [] for _ in range(num_experiments):     p_val, effect_size = simulate_experiment()     p_values.append(p_val)     effect_sizes.append(effect_size) # 只包括統計上顯著的結果過濾 significant_effect_sizes = [effect_sizes[i] for i ...

A New Look at p-Values for Randomized Clinical Trials 23

#  https://t.co/iIu46aY3mZ # http://www.stat.columbia.edu/~gelman/research/published/pval_RCTs_rev3.pdf # Evaluating a shrinkage estimator for the treatment effect in clinical trials 23 # https://onlinelibrary.wiley.com/doi/10.1002/sim.9992 # b= β+N(0, s²), SNR= β/σ, z=SNR+N(0, 1), β true effect # We compute the shrinkage estimator β̂ and its SD σ from the unbiased estimator b and its SD s with the following R code: shrinkage <- function (b, s) { z <- b/s p<- c(0.32, 0.31, 0.3, 0.07) # from Table 3 sigma <- c(0.61, 1.42, 2.16, 5.64) sigma2 <- sigma^2 q <-p*dnorm (z, 0, sqrt (sigma2+1)) q <-q/sum (q) # conditional mixing probs pm <-b*sigma2/(sigma2+1) # conditional means pv <- s^2*sigma2/(sigma2+1) # conditional variances data.frame(q, pm, pv) } # For example, if we observe b = 0.4 and s = 0.3 then the usual 95% confidence interval is -0.2 to 1.0. We can compute the shrinkage estimator as follows: shrink <- shrinkage (b=0.4, s=0.3) betahat <- sum (sh...