RSK World - Statsmodels Statistical Modeling - Project Files | RSK World

bayesian_statistics.py performance_benchmarking.py

bayesian_statistics.py

"""
Bayesian Statistical Analysis

Author: RSK World
Website: https://rskworld.in
Email: help@rskworld.in
Phone: +91 93305 39277
"""

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
import warnings
warnings.filterwarnings('ignore')


class BayesianAnalysis:
    """
    Bayesian Statistical Analysis Tools
    
    Author: RSK World
    Website: https://rskworld.in
    Email: help@rskworld.in
    Phone: +91 93305 39277
    """
    
    @staticmethod
    def bayesian_ttest(sample1, sample2, prior_mean=0, prior_std=1):
        """
        Bayesian t-test
        
        Parameters:
        -----------
        sample1 : array-like
            First sample
        sample2 : array-like
            Second sample
        prior_mean : float
            Prior mean
        prior_std : float
            Prior standard deviation
        """
        n1, n2 = len(sample1), len(sample2)
        mean1, mean2 = np.mean(sample1), np.mean(sample2)
        var1, var2 = np.var(sample1, ddof=1), np.var(sample2, ddof=1)
        
        # Pooled variance
        pooled_var = ((n1 - 1) * var1 + (n2 - 1) * var2) / (n1 + n2 - 2)
        pooled_std = np.sqrt(pooled_var)
        
        # Standard error
        se = pooled_std * np.sqrt(1/n1 + 1/n2)
        
        # Posterior distribution parameters
        posterior_mean = (prior_mean / prior_std**2 + (mean1 - mean2) / se**2) / \
                        (1 / prior_std**2 + 1 / se**2)
        posterior_var = 1 / (1 / prior_std**2 + 1 / se**2)
        posterior_std = np.sqrt(posterior_var)
        
        # Credible interval (95%)
        ci_lower = posterior_mean - 1.96 * posterior_std
        ci_upper = posterior_mean + 1.96 * posterior_std
        
        print("Bayesian t-test Results:")
        print("=" * 70)
        print(f"Posterior Mean: {posterior_mean:.4f}")
        print(f"Posterior Std: {posterior_std:.4f}")
        print(f"95% Credible Interval: [{ci_lower:.4f}, {ci_upper:.4f}]")
        
        # Probability that difference > 0
        prob_positive = 1 - stats.norm.cdf(0, posterior_mean, posterior_std)
        print(f"Probability (difference > 0): {prob_positive:.4f}")
        
        return {
            'posterior_mean': posterior_mean,
            'posterior_std': posterior_std,
            'ci_lower': ci_lower,
            'ci_upper': ci_upper,
            'prob_positive': prob_positive
        }
    
    @staticmethod
    def bayesian_linear_regression(X, y, prior_precision=1.0):
        """
        Bayesian Linear Regression
        
        Parameters:
        -----------
        X : array-like
            Independent variables
        y : array-like
            Dependent variable
        prior_precision : float
            Prior precision (inverse variance)
        """
        from statsmodels.api import add_constant, OLS
        
        X_with_const = add_constant(X)
        n, p = X_with_const.shape
        
        # OLS estimates
        ols_model = OLS(y, X_with_const).fit()
        beta_ols = ols_model.params
        sigma_sq = ols_model.mse_resid
        
        # Prior covariance
        prior_cov = np.eye(p) / prior_precision
        
        # Posterior covariance
        XtX = X_with_const.T @ X_with_const
        posterior_cov = np.linalg.inv(prior_precision * np.eye(p) + XtX / sigma_sq)
        
        # Posterior mean
        posterior_mean = posterior_cov @ (XtX @ beta_ols / sigma_sq)
        
        # Credible intervals
        posterior_std = np.sqrt(np.diag(posterior_cov))
        ci_lower = posterior_mean - 1.96 * posterior_std
        ci_upper = posterior_mean + 1.96 * posterior_std
        
        print("Bayesian Linear Regression Results:")
        print("=" * 70)
        results_df = pd.DataFrame({
            'Posterior Mean': posterior_mean,
            'Posterior Std': posterior_std,
            'CI Lower': ci_lower,
            'CI Upper': ci_upper
        }, index=[f'beta_{i}' for i in range(p)])
        print(results_df)
        
        return {
            'posterior_mean': posterior_mean,
            'posterior_cov': posterior_cov,
            'posterior_std': posterior_std,
            'ci_lower': ci_lower,
            'ci_upper': ci_upper
        }
    
    @staticmethod
    def plot_posterior_distribution(mean, std, true_value=None, title="Posterior Distribution"):
        """Plot posterior distribution"""
        x = np.linspace(mean - 4*std, mean + 4*std, 1000)
        y = stats.norm.pdf(x, mean, std)
        
        plt.figure(figsize=(10, 6))
        plt.plot(x, y, 'b-', linewidth=2, label='Posterior')
        plt.axvline(mean, color='r', linestyle='--', label=f'Mean: {mean:.4f}')
        
        if true_value is not None:
            plt.axvline(true_value, color='g', linestyle='--', label=f'True: {true_value:.4f}')
        
        # Credible interval
        ci_lower = mean - 1.96 * std
        ci_upper = mean + 1.96 * std
        plt.axvspan(ci_lower, ci_upper, alpha=0.2, color='blue', label='95% CI')
        
        plt.xlabel('Parameter Value')
        plt.ylabel('Density')
        plt.title(title)
        plt.legend()
        plt.grid(True, alpha=0.3)
        plt.tight_layout()
        plt.show()


class BayesianModelComparison:
    """
    Bayesian Model Comparison
    
    Author: RSK World
    Website: https://rskworld.in
    Email: help@rskworld.in
    Phone: +91 93305 39277
    """
    
    @staticmethod
    def bayes_factor(model1_loglik, model2_loglik):
        """
        Calculate Bayes Factor
        
        Parameters:
        -----------
        model1_loglik : float
            Log-likelihood of model 1
        model2_loglik : float
            Log-likelihood of model 2
        """
        bf = np.exp(model1_loglik - model2_loglik)
        
        print("Bayes Factor:")
        print("=" * 70)
        print(f"BF (Model 1 / Model 2): {bf:.4f}")
        
        if bf > 1:
            print(f"Model 1 is {bf:.2f}x more likely than Model 2")
        else:
            print(f"Model 2 is {1/bf:.2f}x more likely than Model 1")
        
        # Interpretation
        if bf > 100:
            strength = "Very strong evidence"
        elif bf > 10:
            strength = "Strong evidence"
        elif bf > 3:
            strength = "Moderate evidence"
        else:
            strength = "Weak evidence"
        
        print(f"Evidence strength: {strength}")
        
        return bf
    
    @staticmethod
    def bayesian_information_criterion(loglik, n_params, n_obs):
        """
        Calculate Bayesian Information Criterion (BIC)
        
        Parameters:
        -----------
        loglik : float
            Log-likelihood
        n_params : int
            Number of parameters
        n_obs : int
            Number of observations
        """
        bic = -2 * loglik + n_params * np.log(n_obs)
        return bic


if __name__ == "__main__":
    # Example usage
    print("Bayesian Statistics Example")
    print("=" * 70)
    
    # Generate sample data
    np.random.seed(42)
    sample1 = np.random.normal(100, 15, 30)
    sample2 = np.random.normal(105, 15, 30)
    
    # Bayesian t-test
    result = BayesianAnalysis.bayesian_ttest(sample1, sample2)
    
    # Plot posterior
    BayesianAnalysis.plot_posterior_distribution(
        result['posterior_mean'], 
        result['posterior_std'],
        title="Posterior Distribution of Mean Difference"
    )

252 lines•7.5 KB

python

performance_benchmarking.py

Raw Download

"""
Performance Benchmarking and Profiling

Author: RSK World
Website: https://rskworld.in
Email: help@rskworld.in
Phone: +91 93305 39277
"""

import time
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from functools import wraps
import warnings
warnings.filterwarnings('ignore')


class PerformanceBenchmark:
    """
    Performance Benchmarking Tools
    
    Author: RSK World
    Website: https://rskworld.in
    Email: help@rskworld.in
    Phone: +91 93305 39277
    """
    
    def __init__(self):
        self.benchmarks = []
    
    def benchmark_function(self, func, *args, **kwargs):
        """
        Benchmark a function execution
        
        Parameters:
        -----------
        func : callable
            Function to benchmark
        *args, **kwargs
            Arguments to pass to function
        """
        start_time = time.time()
        start_memory = self._get_memory_usage()
        
        result = func(*args, **kwargs)
        
        end_time = time.time()
        end_memory = self._get_memory_usage()
        
        execution_time = end_time - start_time
        memory_used = end_memory - start_memory
        
        benchmark = {
            'function': func.__name__,
            'execution_time': execution_time,
            'memory_used': memory_used,
            'timestamp': time.time()
        }
        
        self.benchmarks.append(benchmark)
        
        print(f"Benchmark: {func.__name__}")
        print(f"Execution Time: {execution_time:.4f} seconds")
        print(f"Memory Used: {memory_used:.2f} MB")
        
        return result, benchmark
    
    def compare_models(self, models_dict, X, y):
        """
        Compare performance of multiple models
        
        Parameters:
        -----------
        models_dict : dict
            Dictionary of model names and model objects
        X : array-like
            Independent variables
        y : array-like
            Dependent variable
        """
        results = []
        
        for name, model in models_dict.items():
            print(f"\nBenchmarking: {name}")
            print("-" * 70)
            
            # Training time
            start = time.time()
            if hasattr(model, 'fit'):
                model.fit(X, y)
            training_time = time.time() - start
            
            # Prediction time
            start = time.time()
            if hasattr(model, 'predict'):
                predictions = model.predict(X)
            prediction_time = time.time() - start
            
            # Model metrics
            metrics = {}
            if hasattr(model, 'results'):
                metrics['aic'] = model.results.aic if hasattr(model.results, 'aic') else None
                metrics['bic'] = model.results.bic if hasattr(model.results, 'bic') else None
                metrics['r_squared'] = model.results.rsquared if hasattr(model.results, 'rsquared') else None
            
            results.append({
                'model': name,
                'training_time': training_time,
                'prediction_time': prediction_time,
                'total_time': training_time + prediction_time,
                **metrics
            })
        
        results_df = pd.DataFrame(results)
        
        print("\n" + "=" * 70)
        print("MODEL PERFORMANCE COMPARISON")
        print("=" * 70)
        print(results_df.to_string(index=False))
        
        return results_df
    
    def plot_benchmark_comparison(self, metric='execution_time'):
        """Plot benchmark comparison"""
        if not self.benchmarks:
            print("No benchmarks to plot")
            return
        
        df = pd.DataFrame(self.benchmarks)
        
        plt.figure(figsize=(10, 6))
        plt.bar(df['function'], df[metric])
        plt.xlabel('Function')
        plt.ylabel(metric.replace('_', ' ').title())
        plt.title('Performance Benchmark Comparison')
        plt.xticks(rotation=45, ha='right')
        plt.grid(True, alpha=0.3, axis='y')
        plt.tight_layout()
        plt.show()
    
    def _get_memory_usage(self):
        """Get current memory usage in MB"""
        try:
            import psutil
            process = psutil.Process()
            return process.memory_info().rss / 1024 / 1024  # Convert to MB
        except ImportError:
            return 0
    
    def timeit(self, n_iterations=10):
        """
        Decorator for timing function execution
        
        Parameters:
        -----------
        n_iterations : int
            Number of iterations to average
        """
        def decorator(func):
            @wraps(func)
            def wrapper(*args, **kwargs):
                times = []
                for _ in range(n_iterations):
                    start = time.time()
                    result = func(*args, **kwargs)
                    times.append(time.time() - start)
                
                avg_time = np.mean(times)
                std_time = np.std(times)
                
                print(f"{func.__name__}: {avg_time:.4f} ± {std_time:.4f} seconds (avg of {n_iterations} runs)")
                
                return result
            return wrapper
        return decorator


def benchmark_model_complexity(X, y, model_class, complexity_params):
    """
    Benchmark model performance across different complexities
    
    Parameters:
    -----------
    X : array-like
        Independent variables
    y : array-like
        Dependent variable
    model_class : class
        Model class to instantiate
    complexity_params : dict
        Dictionary of complexity parameter names and values
    """
    results = []
    
    for param_name, param_values in complexity_params.items():
        for param_value in param_values:
            model = model_class()
            setattr(model, param_name, param_value)
            
            start = time.time()
            model.fit(X, y)
            training_time = time.time() - start
            
            start = time.time()
            predictions = model.predict(X)
            prediction_time = time.time() - start
            
            results.append({
                param_name: param_value,
                'training_time': training_time,
                'prediction_time': prediction_time,
                'total_time': training_time + prediction_time
            })
    
    return pd.DataFrame(results)


if __name__ == "__main__":
    # Example usage
    print("Performance Benchmarking Example")
    print("=" * 70)
    
    from regression_analysis import LinearRegressionModel
    import numpy as np
    
    # Generate sample data
    np.random.seed(42)
    X = np.random.randn(1000, 10)
    y = 2 + np.sum(X[:, :3], axis=1) + np.random.randn(1000) * 0.5
    
    # Benchmark
    benchmark = PerformanceBenchmark()
    
    def fit_model():
        model = LinearRegressionModel()
        model.fit(X, y)
        return model
    
    result, bench = benchmark.benchmark_function(fit_model)
    
    # Compare models
    models = {
        'Simple Model': LinearRegressionModel()
    }
    
    comparison = benchmark.compare_models(models, X, y)

249 lines•7.3 KB

python

Theme Settings

Color Scheme

Display Options

Font Size

Theme Settings

Color Scheme

Display Options

Font Size

🚀 Support RSK World