#!/usr/bin/env python3
"""
Chapter 12: Python for Next Gen Biologist
Script: ch12_ml_cancer_classifier.py
Description: Machine Learning for Biology:
             Feature scaling, Principal Component Analysis (PCA),
             and Random Forest Classifier for cancer subtype classification.
Usage:
    python ch12_ml_cancer_classifier.py
"""

import numpy as np

def run_ml_pipeline():
    try:
        from sklearn.ensemble import RandomForestClassifier
        from sklearn.decomposition import PCA
        from sklearn.model_selection import StratifiedKFold, cross_val_score
        from sklearn.preprocessing import StandardScaler
        
        print("=== Chapter 12: Machine Learning for Biology (Scikit-Learn) ===")
        np.random.seed(42)
        
        # Simulate 60 patient expression profiles (30 Luminal A, 30 Basal-like) across 200 biomarker genes
        n_samples = 60
        n_features = 200
        
        X = np.random.normal(loc=5.0, scale=1.0, size=(n_samples, n_features))
        y = np.array([0] * 30 + [1] * 30)  # 0 = Luminal A, 1 = Basal-like
        
        # Introduce distinct biomarker signature in first 10 genes
        X[:30, :10] += 2.0   # High in Luminal A
        X[30:, 10:20] += 2.0 # High in Basal-like
        
        # 1. Feature Standardization
        scaler = StandardScaler()
        X_scaled = scaler.fit_transform(X)
        
        # 2. PCA Dimensionality Reduction
        pca = PCA(n_components=2)
        X_pca = pca.fit_transform(X_scaled)
        var_exp = pca.explained_variance_ratio_ * 100
        print(f"\n[1] PCA Decomposition:")
        print(f"    - PC1 Variance Explained: {var_exp[0]:.2f}%")
        print(f"    - PC2 Variance Explained: {var_exp[1]:.2f}%")
        
        # 3. Random Forest Classification with 5-Fold Stratified Cross-Validation
        rf = RandomForestClassifier(n_estimators=50, random_state=42)
        cv = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
        scores = cross_val_score(rf, X_scaled, y, cv=cv, scoring='accuracy')
        
        print(f"\n[2] Random Forest 5-Fold Cross-Validation Accuracy:")
        print(f"    - Mean Accuracy : {np.mean(scores)*100:.2f}% (+/- {np.std(scores)*100:.2f}%)")
        print(f"    - Fold Scores   : {[round(s, 3) for s in scores]}")
        
        # 4. Feature Importance for Biomarker Discovery
        rf.fit(X_scaled, y)
        importances = rf.feature_importances_
        top_genes = np.argsort(importances)[::-1][:5]
        print(f"\n[3] Top 5 Predictive Biomarker Genes:")
        for rank, gene_idx in enumerate(top_genes, 1):
            print(f"    Rank {rank}: Gene_{gene_idx:03d} (Importance: {importances[gene_idx]:.4f})")
            
    except ImportError:
        print("[NOTICE] 'scikit-learn' is not installed in current environment.")
        print("Run: pip install scikit-learn numpy pandas")

def main():
    run_ml_pipeline()

if __name__ == "__main__":
    main()