example_statistical_analysis.py

"""
Example: Using Statistical Analysis Functions

This demonstrates how to use the new statistical analysis features
added to tests_paper_analysis.py for paper writing.
"""

import json
import sys
from pathlib import Path

sys.path.insert(0, str(Path(__file__).parent / 'src'))

from tests_paper_analysis import (
    analyze_results,
    compute_statistical_tests,
    compute_cohens_d,
    compute_confidence_interval
)
import numpy as np


def example_full_analysis():
    """Example: Run full analysis pipeline"""
    print("\n" + "=" * 80)
    print("EXAMPLE 1: Full Analysis Pipeline")
    print("=" * 80)
    
    # Run complete analysis
    results = analyze_results(
        test_results_json='test_results/test_results.json',
        output_dir='/tmp/paper_figures/',
        verbose=True
    )
    
    print("\nGenerated files in /tmp/paper_figures/")
    print("- statistical_report.txt (includes p-values, effect sizes, CIs)")


def example_statistical_tests():
    """Example: Use statistical tests directly"""
    print("\n" + "=" * 80)
    print("EXAMPLE 2: Statistical Tests on Test Results")
    print("=" * 80)
    
    # Load test results
    with open('test_results/test_results.json', 'r') as f:
        results_dict = json.load(f)
    
    # Compute statistical tests
    stats = compute_statistical_tests(results_dict)
    
    print("\nStatistical Test Results:")
    print("-" * 80)
    print(f"Control Mean (no residue):     {stats['control_mean']:.4f}")
    print(f"Treatment Mean (with residue): {stats['treatment_mean']:.4f}")
    
    if not np.isnan(stats['control_std']):
        print(f"Control Std Dev:               {stats['control_std']:.4f}")
    if not np.isnan(stats['treatment_std']):
        print(f"Treatment Std Dev:             {stats['treatment_std']:.4f}")
    
    if stats.get('t_stat') is not None:
        print(f"\nT-statistic:                   {stats['t_stat']:.4f}")
        print(f"P-value:                       {stats['p_value']:.4f}")
        
        # Interpret significance
        p = stats['p_value']
        if p < 0.001:
            print("Significance:                  *** (p < 0.001, highly significant)")
        elif p < 0.01:
            print("Significance:                  ** (p < 0.01, very significant)")
        elif p < 0.05:
            print("Significance:                  * (p < 0.05, significant)")
        else:
            print("Significance:                  n.s. (p >= 0.05, not significant)")
    
    if stats.get('effect_size') is not None:
        print(f"\nCohen's d (effect size):       {stats['effect_size']:.4f}")
        print(f"Interpretation:                {stats['interpretation']}")
        
        # Provide context
        d = abs(stats['effect_size'])
        if d < 0.2:
            print("→ Negligible practical difference between groups")
        elif d < 0.5:
            print("→ Small but noticeable difference")
        elif d < 0.8:
            print("→ Medium, meaningful difference")
        else:
            print("→ Large, highly meaningful difference")
    
    if 'warning' in stats:
        print(f"\nNote: {stats['warning']}")


def example_custom_data():
    """Example: Compute statistics on custom data"""
    print("\n" + "=" * 80)
    print("EXAMPLE 3: Custom Data Analysis")
    print("=" * 80)
    
    # Example: Harm rates from different experimental runs
    control_harm_rates = np.array([0.45, 0.42, 0.48, 0.44, 0.46])
    treatment_harm_rates = np.array([0.25, 0.22, 0.28, 0.24, 0.26])
    
    # Compute effect size
    d = compute_cohens_d(control_harm_rates, treatment_harm_rates)
    
    # Compute confidence intervals
    control_ci = compute_confidence_interval(control_harm_rates)
    treatment_ci = compute_confidence_interval(treatment_harm_rates)
    
    print("\nHarm Rate Analysis:")
    print("-" * 80)
    print(f"Control Mean:    {np.mean(control_harm_rates):.4f}")
    print(f"Control 95% CI:  [{control_ci[0]:.4f}, {control_ci[1]:.4f}]")
    print(f"\nTreatment Mean:  {np.mean(treatment_harm_rates):.4f}")
    print(f"Treatment 95% CI: [{treatment_ci[0]:.4f}, {treatment_ci[1]:.4f}]")
    print(f"\nCohen's d:       {d:.4f}")
    
    print("\n" + "=" * 80)
    print("FOR PAPER (APA Style):")
    print("=" * 80)
    
    # Compute means once
    control_mean = np.mean(control_harm_rates)
    treatment_mean = np.mean(treatment_harm_rates)
    
    # Determine direction of effect
    if treatment_mean < control_mean:
        effect_desc = "reduced harm rates"
    else:
        effect_desc = "changed harm rates"
    
    print(f"The residue field {effect_desc} from M = {control_mean:.2f}")
    print(f"(95% CI [{control_ci[0]:.2f}, {control_ci[1]:.2f}]) to M = {treatment_mean:.2f}")
    print(f"(95% CI [{treatment_ci[0]:.2f}, {treatment_ci[1]:.2f}]), a large effect")
    print(f"(Cohen's d = {d:.2f}).")


def example_paper_reporting():
    """Example: Extract statistics for paper writing"""
    print("\n" + "=" * 80)
    print("EXAMPLE 4: Paper Reporting Template")
    print("=" * 80)
    
    # Load results
    with open('test_results/test_results.json', 'r') as f:
        results_dict = json.load(f)
    
    stats = compute_statistical_tests(results_dict)
    
    print("\nResults Section Template:")
    print("-" * 80)
    print(f"""
We compared harmful action rates between control (no residue field) and 
treatment (with residue field) conditions. The treatment group showed 
{'higher' if stats['treatment_mean'] > stats['control_mean'] else 'lower'} 
harm rates (M = {stats['treatment_mean']:.2f}) compared to control 
(M = {stats['control_mean']:.2f}).""")
    
    if stats.get('t_stat') is not None:
        print(f"""
This difference was {'statistically significant' if stats['p_value'] < 0.05 else 'not statistically significant'}
(t = {stats['t_stat']:.2f}, p {'< 0.001' if stats['p_value'] < 0.001 else f'= {stats["p_value"]:.3f}'}).""")
    
    if stats.get('effect_size') is not None:
        print(f"""
The effect size was {stats['interpretation']} (Cohen's d = {stats['effect_size']:.2f}),
indicating a {'negligible' if abs(stats['effect_size']) < 0.2 else 'meaningful'} 
practical difference between conditions.
    """)


if __name__ == '__main__':
    print("\n" + "=" * 80)
    print("STATISTICAL ANALYSIS EXAMPLES")
    print("=" * 80)
    
    try:
        # Only run statistical tests example to avoid file I/O issues
        example_statistical_tests()
        example_custom_data()
        example_paper_reporting()
        
        print("\n" + "=" * 80)
        print("✓ All examples completed successfully")
        print("=" * 80)
        print("\nFor full analysis with plots, run:")
        print("  python tests_paper_analysis.py --input test_results/test_results.json --output figures/")
        print("\n")
        
    except Exception as e:
        print(f"\n✗ Error: {e}\n")
        import traceback
        traceback.print_exc()
        sys.exit(1)