Winner’s curse for gaussian dist

# Simulate winner’s curse of choosing only those with p of less than 0.05 of gaussian dist with t tests in clinical trials using python

import numpy as np

from scipy import stats

import matplotlib.pyplot as plt

# Setting the seed for reproducibility

np.random.seed(42)

# Parameters

true_effect = 0.5  # True mean difference between treatment and placebo

std_dev = 3  # Standard deviation for both groups

num_trials = 1000  # Number of trials to simulate

sample_size = 300  # Number of participants in each group per trial

# Store results

p_values = []

effect_sizes = []

for _ in range(num_trials):

    # Generate data for treatment and placebo groups

    treatment_group = np.random.normal(true_effect, std_dev, sample_size)

    placebo_group = np.random.normal(0, std_dev, sample_size)

    # Calculate the effect size (mean difference)

    effect_size = np.mean(treatment_group) - np.mean(placebo_group)

    effect_sizes.append(effect_size)

    # Perform a t-test

    t_stat, p_val = stats.ttest_ind(treatment_group, placebo_group)

    p_values.append(p_val)

# Filter results where p-value < 0.05

significant_effect_sizes = [effect_sizes[i] for i in range(num_trials) if p_values[i] < 0.05]

# Analysis

print(f"True effect: {true_effect}")

print(f"Mean effect size across all trials: {np.mean(effect_sizes)}")

print(f"Mean effect size of significant trials (p < 0.05): {np.mean(significant_effect_sizes) if significant_effect_sizes else 'No significant results'}")

# Plotting

plt.figure(figsize=(12, 6))

plt.hist(effect_sizes, bins=20, color='blue', alpha=0.7, label='All Trials')

plt.hist(significant_effect_sizes, bins=20, color='red', alpha=0.7, label=f'Significant Trials (p < 0.05), n={len(significant_effect_sizes)}, mean={np.mean(significant_effect_sizes):.2f}')

plt.axvline(x=true_effect, color='green', linestyle='--', label='True Effect')

plt.axvline(x=np.mean(significant_effect_sizes) if significant_effect_sizes else 0, color='red', linestyle='--', label=f'Mean of Significant Trials')

plt.title('Effect Sizes of Trials')

plt.xlabel('Effect Size')

plt.ylabel('Frequency')

plt.legend()

plt.show()