{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.11.11","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":96164,"databundleVersionId":11418275,"sourceType":"competition"}],"dockerImageVersionId":31041,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2025-06-17T05:23:08.981627Z","iopub.execute_input":"2025-06-17T05:23:08.981829Z","iopub.status.idle":"2025-06-17T05:23:10.823089Z","shell.execute_reply.started":"2025-06-17T05:23:08.981812Z","shell.execute_reply":"2025-06-17T05:23:10.822294Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#!/usr/bin/env python3\n\"\"\"\nComprehensive DRW Crypto Dimensionality Reduction Analysis - Enhanced Sample Size\n===============================================================================\n\nComplete evaluation of 25+ dimensionality reduction techniques with robust preprocessing\nUsing larger sample sizes for more comprehensive analysis\n\"\"\"\n\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport warnings\nimport os\nwarnings.filterwarnings('ignore')\n\n# Core ML libraries\nfrom sklearn.preprocessing import StandardScaler, RobustScaler, MinMaxScaler, PowerTransformer, QuantileTransformer\nfrom sklearn.decomposition import PCA, KernelPCA, FactorAnalysis, FastICA, NMF, DictionaryLearning\nfrom sklearn.decomposition import TruncatedSVD, SparsePCA, MiniBatchSparsePCA, IncrementalPCA\nfrom sklearn.decomposition import LatentDirichletAllocation, MiniBatchDictionaryLearning\nfrom sklearn.manifold import TSNE, Isomap, LocallyLinearEmbedding, MDS, SpectralEmbedding\nfrom sklearn.discriminant_analysis import LinearDiscriminantAnalysis\nfrom sklearn.random_projection import GaussianRandomProjection, SparseRandomProjection\nfrom sklearn.feature_selection import SelectKBest, f_regression, mutual_info_regression, VarianceThreshold\nfrom sklearn.metrics import mean_squared_error, r2_score\nfrom sklearn.ensemble import RandomForestRegressor, GradientBoostingRegressor\nfrom sklearn.linear_model import ElasticNet, Ridge, Lasso\nfrom sklearn.model_selection import train_test_split, cross_val_score\nfrom sklearn.neighbors import NeighborhoodComponentsAnalysis\nfrom sklearn.cross_decomposition import PLSRegression, CCA\nfrom sklearn.pipeline import Pipeline\nfrom scipy import stats\nfrom scipy.stats import pearsonr, mstats\nimport time\nfrom datetime import datetime\n\n# Advanced libraries\ntry:\n    import umap\n    UMAP_AVAILABLE = True\nexcept ImportError:\n    UMAP_AVAILABLE = False\n    print(\"UMAP not available - UMAP methods will be skipped\")\n\ntry:\n    import tensorflow as tf\n    from tensorflow import keras\n    from tensorflow.keras import layers, Model, regularizers\n    tf.random.set_seed(42)\n    TENSORFLOW_AVAILABLE = True\nexcept ImportError:\n    TENSORFLOW_AVAILABLE = False\n    print(\"TensorFlow not available - deep learning methods will be skipped\")\n\ntry:\n    from minisom import MiniSom\n    MINISOM_AVAILABLE = True\nexcept ImportError:\n    MINISOM_AVAILABLE = False\n\n# Set random seed\nnp.random.seed(42)\n\n# Define important features from your analysis\nIMPORTANT_FEATURES = [\n    \"X863\", \"X856\", \"X344\", \"X598\", \"X862\", \"X385\", \"X852\", \"X603\", \"X860\", \"X674\",\n    \"X415\", \"X345\", \"X137\", \"X855\", \"X174\", \"X302\", \"X178\", \"X532\", \"X168\", \"X612\",\n    \"bid_qty\", \"ask_qty\", \"buy_qty\", \"sell_qty\", \"volume\", \"X888\", \"X421\", \"X333\"\n]\n\n# Check if we're in Kaggle environment\nKAGGLE_INPUT_PATH = '/kaggle/input/drw-crypto-market-prediction'\nif os.path.exists(KAGGLE_INPUT_PATH):\n    DATA_PATH = KAGGLE_INPUT_PATH\nelse:\n    DATA_PATH = '.'\n\n# INCREASED SAMPLE SIZES\nDEFAULT_SAMPLE_SIZE = 25000  # Increased from 10000\nSYNTHETIC_SAMPLE_SIZE = 25000  # Increased from 5000\nPREPROCESSING_EVAL_SIZE = 25000  # Increased from 5000\n\n\nclass RobustDataPreprocessor:\n    \"\"\"Enhanced preprocessor with multiple strategies for handling infinite values\"\"\"\n    \n    def __init__(self):\n        self.preprocessing_strategy = None\n        self.feature_scalers = {}\n        self.outlier_bounds = {}\n        self.variance_selector = None\n        self.winsorize_limits = {}\n        self.feature_stats = {}\n        \n    def detect_inf_features(self, X, feature_names):\n        \"\"\"Detect features with infinite values\"\"\"\n        inf_features = []\n        for i, name in enumerate(feature_names):\n            if np.any(np.isinf(X[:, i])):\n                n_pos_inf = np.sum(np.isposinf(X[:, i]))\n                n_neg_inf = np.sum(np.isneginf(X[:, i]))\n                inf_features.append({\n                    'index': i,\n                    'name': name,\n                    'pos_inf': n_pos_inf,\n                    'neg_inf': n_neg_inf,\n                    'total_inf': n_pos_inf + n_neg_inf,\n                    'pct_inf': (n_pos_inf + n_neg_inf) / len(X) * 100\n                })\n        return inf_features\n    \n    def winsorize_features(self, X, limits=(0.01, 0.99)):\n        \"\"\"Apply winsorization to handle extreme values\"\"\"\n        X_winsorized = X.copy()\n        \n        for i in range(X.shape[1]):\n            col = X[:, i]\n            \n            # Handle infinite values first\n            finite_mask = np.isfinite(col)\n            if np.any(finite_mask):\n                finite_vals = col[finite_mask]\n                \n                # Calculate winsorization limits on finite values\n                lower_limit = np.percentile(finite_vals, limits[0] * 100)\n                upper_limit = np.percentile(finite_vals, limits[1] * 100)\n                \n                # Replace infinities\n                col[np.isneginf(col)] = lower_limit\n                col[np.isposinf(col)] = upper_limit\n                \n                # Winsorize the rest\n                col = np.clip(col, lower_limit, upper_limit)\n                \n                self.winsorize_limits[i] = (lower_limit, upper_limit)\n            else:\n                # All values are infinite - replace with 0\n                col[:] = 0\n                self.winsorize_limits[i] = (0, 0)\n            \n            X_winsorized[:, i] = col\n            \n        return X_winsorized\n    \n    def quartile_capping(self, X, iqr_multiplier=3.0):\n        \"\"\"Cap outliers using IQR method\"\"\"\n        X_capped = X.copy()\n        \n        for i in range(X.shape[1]):\n            col = X[:, i]\n            \n            # Work with finite values\n            finite_mask = np.isfinite(col)\n            if np.any(finite_mask):\n                finite_vals = col[finite_mask]\n                \n                Q1 = np.percentile(finite_vals, 25)\n                Q3 = np.percentile(finite_vals, 75)\n                IQR = Q3 - Q1\n                \n                lower_bound = Q1 - iqr_multiplier * IQR\n                upper_bound = Q3 + iqr_multiplier * IQR\n                \n                # Replace infinities\n                col[np.isneginf(col)] = lower_bound\n                col[np.isposinf(col)] = upper_bound\n                \n                # Cap the rest\n                col = np.clip(col, lower_bound, upper_bound)\n                \n                self.outlier_bounds[i] = (lower_bound, upper_bound)\n            else:\n                col[:] = 0\n                self.outlier_bounds[i] = (0, 0)\n            \n            X_capped[:, i] = col\n            \n        return X_capped\n    \n    def percentile_clipping(self, X, lower_pct=0.1, upper_pct=99.9):\n        \"\"\"Clip values to percentile range\"\"\"\n        X_clipped = X.copy()\n        \n        for i in range(X.shape[1]):\n            col = X[:, i]\n            \n            # Work with finite values\n            finite_mask = np.isfinite(col)\n            if np.any(finite_mask):\n                finite_vals = col[finite_mask]\n                \n                lower_val = np.percentile(finite_vals, lower_pct)\n                upper_val = np.percentile(finite_vals, upper_pct)\n                \n                # Replace infinities\n                col[np.isneginf(col)] = lower_val\n                col[np.isposinf(col)] = upper_val\n                \n                # Clip the rest\n                col = np.clip(col, lower_val, upper_val)\n            else:\n                col[:] = 0\n            \n            X_clipped[:, i] = col\n            \n        return X_clipped\n    \n    def select_best_preprocessing(self, X, y, feature_names):\n        \"\"\"Test multiple preprocessing strategies and select the best\"\"\"\n        \n        strategies = [\n            {\n                'name': 'winsorize_01_99',\n                'method': lambda X: self.winsorize_features(X, limits=(0.01, 0.99))\n            },\n            {\n                'name': 'quartile_cap_3',\n                'method': lambda X: self.quartile_capping(X, iqr_multiplier=3.0)\n            },\n            {\n                'name': 'percentile_clip_0.1_99.9',\n                'method': lambda X: self.percentile_clipping(X, lower_pct=0.1, upper_pct=99.9)\n            }\n        ]\n        \n        # Detect features with infinite values\n        inf_features = self.detect_inf_features(X, feature_names)\n        if inf_features:\n            print(f\"\\nFound {len(inf_features)} features with infinite values\")\n            # Show summary of top features with infinities\n            inf_summary = sorted(inf_features, key=lambda x: x['pct_inf'], reverse=True)[:5]\n            for feat in inf_summary:\n                print(f\"  {feat['name']}: {feat['pct_inf']:.2f}% infinite values\")\n        \n        best_score = -np.inf\n        best_strategy = strategies[0]\n        \n        print(\"\\nTesting preprocessing strategies...\")\n        for strategy in strategies:\n            try:\n                # Apply preprocessing\n                X_processed = strategy['method'](X.copy())\n                \n                # Check if any infinities remain\n                n_inf_remaining = np.sum(np.isinf(X_processed))\n                \n                if n_inf_remaining > 0:\n                    print(f\"  {strategy['name']}: Still has {n_inf_remaining} infinite values - skipping\")\n                    continue\n                \n                # Scale data\n                scaler = RobustScaler()\n                X_scaled = scaler.fit_transform(X_processed)\n                \n                # Quick evaluation with simple model - use larger sample\n                rf = RandomForestRegressor(n_estimators=50, max_depth=8, random_state=42, n_jobs=-1)\n                \n                # Use larger subset for evaluation\n                n_eval = min(PREPROCESSING_EVAL_SIZE, len(X_scaled))\n                print(f\"  Evaluating {strategy['name']} on {n_eval} samples...\")\n                \n                scores = cross_val_score(rf, X_scaled[:n_eval], y[:n_eval], cv=3, scoring='r2')\n                score = np.mean(scores)\n                \n                print(f\"  {strategy['name']}: R² = {score:.4f} (std: {np.std(scores):.4f})\")\n                \n                if score > best_score:\n                    best_score = score\n                    best_strategy = strategy\n                    \n            except Exception as e:\n                print(f\"  {strategy['name']}: Failed - {str(e)}\")\n        \n        print(f\"\\nBest preprocessing strategy: {best_strategy['name']} (R² = {best_score:.4f})\")\n        self.preprocessing_strategy = best_strategy\n        \n        return best_strategy\n    \n    def fit_transform(self, X, y, feature_names):\n        \"\"\"Fit preprocessor and transform data\"\"\"\n        \n        # Select best preprocessing strategy\n        best_strategy = self.select_best_preprocessing(X, y, feature_names)\n        \n        # Apply best preprocessing\n        print(f\"\\nApplying {best_strategy['name']} preprocessing...\")\n        X_processed = best_strategy['method'](X.copy())\n        \n        # Verify no infinities remain\n        n_inf = np.sum(np.isinf(X_processed))\n        n_nan = np.sum(np.isnan(X_processed))\n        print(f\"After preprocessing: {n_inf} infinite values, {n_nan} NaN values\")\n        \n        # Replace any remaining NaN with 0\n        X_processed = np.nan_to_num(X_processed, nan=0.0)\n        \n        # Remove zero variance features\n        self.variance_selector = VarianceThreshold(threshold=1e-10)\n        X_processed = self.variance_selector.fit_transform(X_processed)\n        feature_mask = self.variance_selector.get_support()\n        feature_names_filtered = [f for f, m in zip(feature_names, feature_mask) if m]\n        \n        print(f\"Removed {len(feature_names) - len(feature_names_filtered)} zero-variance features\")\n        \n        # Scale data\n        scaler = RobustScaler()\n        X_scaled = scaler.fit_transform(X_processed)\n        \n        # Final check\n        print(f\"Final shape: {X_scaled.shape}\")\n        print(f\"Final data range: [{np.min(X_scaled):.3f}, {np.max(X_scaled):.3f}]\")\n        \n        return X_scaled, feature_names_filtered\n\n\nclass ExtendedDimensionalityReducer:\n    \"\"\"Extended dimensionality reduction with 25+ methods\"\"\"\n    \n    def __init__(self, target_dims=50, preserve_important_features=True, important_features=None):\n        self.target_dims = target_dims\n        self.preserve_important_features = preserve_important_features\n        self.important_features = important_features or IMPORTANT_FEATURES\n        self.results = {}\n        self.metrics = {}\n        self.fitted_reducers = {}\n        self.timings = {}\n        \n    def evaluate_reduction(self, X_original, X_reduced, y, method_name):\n        \"\"\"Comprehensive evaluation of dimensionality reduction\"\"\"\n        \n        if len(X_reduced) != len(y):\n            print(f\"Warning: {method_name} changed sample size\")\n            return {}\n        \n        # Check for any remaining infinities or NaN\n        if np.any(np.isinf(X_reduced)) or np.any(np.isnan(X_reduced)):\n            print(f\"Warning: {method_name} produced infinite or NaN values\")\n            X_reduced = np.nan_to_num(X_reduced, nan=0.0, posinf=1e6, neginf=-1e6)\n        \n        metrics = {}\n        \n        # Split data\n        X_orig_train, X_orig_test, X_red_train, X_red_test, y_train, y_test = \\\n            train_test_split(X_original, X_reduced, y, test_size=0.3, random_state=42)\n        \n        # Test with simple models for speed\n        models = {\n            'RandomForest': RandomForestRegressor(n_estimators=100, max_depth=10, random_state=42, n_jobs=-1),\n            'Ridge': Ridge(alpha=1.0, random_state=42)\n        }\n        \n        for model_name, model in models.items():\n            try:\n                # Reduced features performance\n                model.fit(X_red_train, y_train)\n                y_pred_red = model.predict(X_red_test)\n                red_pearson, _ = pearsonr(y_test, y_pred_red)\n                \n                metrics[f'{model_name}_reduced_pearson'] = red_pearson\n                \n            except Exception as e:\n                print(f\"Model {model_name} failed for {method_name}: {e}\")\n        \n        # Feature correlations - analyze more features\n        correlations = []\n        n_features_to_analyze = min(X_reduced.shape[1], 200)  # Increased from 100\n        for i in range(n_features_to_analyze):\n            if not (np.any(np.isinf(X_reduced[:, i])) or np.any(np.isnan(X_reduced[:, i]))):\n                corr, _ = pearsonr(X_reduced[:, i], y)\n                if not np.isnan(corr):\n                    correlations.append(abs(corr))\n        \n        metrics['avg_correlation'] = np.mean(correlations) if correlations else 0\n        metrics['max_correlation'] = np.max(correlations) if correlations else 0\n        metrics['std_correlation'] = np.std(correlations) if correlations else 0\n        metrics['n_features'] = X_reduced.shape[1]\n        \n        return metrics\n    \n    def apply_linear_methods(self, X, y, feature_names):\n        \"\"\"Apply linear dimensionality reduction methods\"\"\"\n        \n        print(\"\\n=== LINEAR METHODS ===\")\n        \n        methods = {}\n        \n        # 1. PCA\n        # methods['PCA'] = PCA(n_components=min(self.target_dims, X.shape[1]-1))\n        \n        # 2. Incremental PCA (good for large datasets)\n        methods['Incremental PCA'] = IncrementalPCA(n_components=min(self.target_dims, X.shape[1]-1), \n                                                    batch_size=256)\n        \n        # 3. Truncated SVD\n        methods['Truncated SVD'] = TruncatedSVD(n_components=min(self.target_dims, X.shape[1]-1))\n        \n        # 4. Factor Analysis\n        methods['Factor Analysis'] = FactorAnalysis(n_components=min(self.target_dims, X.shape[1]//2), \n                                                    random_state=42)\n        \n        # 5. FastICA\n        methods['FastICA'] = FastICA(n_components=min(self.target_dims, 50), \n                                    random_state=42, max_iter=200)\n        \n        # 6. Random Projection\n        methods['Random Projection'] = GaussianRandomProjection(\n            n_components=self.target_dims, random_state=42)\n        \n        # 7. Sparse Random Projection\n        methods['Sparse Random Projection'] = SparseRandomProjection(\n            n_components=self.target_dims, random_state=42)\n        \n        # Apply all methods\n        for name, method in methods.items():\n            try:\n                print(f\"Applying {name}...\")\n                start_time = time.time()\n                \n                X_reduced = method.fit_transform(X)\n                \n                self.timings[name] = time.time() - start_time\n                self.results[name] = X_reduced\n                self.fitted_reducers[name] = method\n                self.metrics[name] = self.evaluate_reduction(X, X_reduced, y, name)\n                \n                print(f\"  {name}: {X_reduced.shape[1]} features, time: {self.timings[name]:.2f}s, \"\n                      f\"avg_corr: {self.metrics[name].get('avg_correlation', 0):.4f}\")\n                    \n            except Exception as e:\n                print(f\"  {name} failed: {e}\")\n    \n    def apply_nonlinear_manifold_methods(self, X, y, feature_names):\n        \"\"\"Apply non-linear manifold learning methods\"\"\"\n        \n        print(\"\\n=== NON-LINEAR MANIFOLD METHODS ===\")\n        \n        methods = {}\n        \n        # 1. Kernel PCA\n        # methods['Kernel PCA (RBF)'] = KernelPCA(n_components=self.target_dims, \n        #                                         kernel='rbf', gamma=0.01,\n        #                                         eigen_solver='randomized')\n        \n        # 2. Kernel PCA with polynomial kernel\n        # methods['Kernel PCA (Poly)'] = KernelPCA(n_components=self.target_dims, \n        #                                          kernel='poly', degree=3,\n        #                                          eigen_solver='randomized')\n        \n        # 3. UMAP (if available)\n        if UMAP_AVAILABLE:\n            n_neighbors = min(30, len(X) // 100)  # Adjusted for larger samples\n            methods['UMAP'] = umap.UMAP(n_components=self.target_dims, \n                                       n_neighbors=n_neighbors, \n                                       min_dist=0.1, random_state=42)\n        \n        # 4. Isomap (limit samples for computational efficiency)\n        if len(X) <= 20000:\n            n_neighbors = min(20, len(X) // 100)\n            methods['Isomap'] = Isomap(n_components=self.target_dims, \n                                      n_neighbors=n_neighbors)\n        \n        # Apply methods\n        for name, method in methods.items():\n            try:\n                print(f\"Applying {name}...\")\n                start_time = time.time()\n                \n                # For computationally expensive methods, use a subset\n                if name in ['Isomap'] and len(X) > 10000:\n                    sample_indices = np.random.choice(len(X), 10000, replace=False)\n                    X_sample = X[sample_indices]\n                    y_sample = y[sample_indices]\n                    X_reduced = method.fit_transform(X_sample)\n                    \n                    # Store reduced results\n                    self.results[name] = X_reduced\n                    self.metrics[name] = self.evaluate_reduction(X_sample, X_reduced, y_sample, name)\n                else:\n                    X_reduced = method.fit_transform(X)\n                    self.results[name] = X_reduced\n                    self.metrics[name] = self.evaluate_reduction(X, X_reduced, y, name)\n                \n                self.timings[name] = time.time() - start_time\n                self.fitted_reducers[name] = method\n                \n                print(f\"  {name}: {X_reduced.shape[1]} features, time: {self.timings[name]:.2f}s, \"\n                      f\"avg_corr: {self.metrics[name].get('avg_correlation', 0):.4f}\")\n                \n            except Exception as e:\n                print(f\"  {name} failed: {e}\")\n    \n    def apply_matrix_factorization_methods(self, X, y, feature_names):\n        \"\"\"Apply matrix factorization methods\"\"\"\n        \n        print(\"\\n=== MATRIX FACTORIZATION METHODS ===\")\n        \n        methods = {}\n        \n        # Ensure non-negative data for NMF\n        X_positive = MinMaxScaler().fit_transform(X)\n        \n        # 1. NMF\n        methods['NMF'] = NMF(n_components=self.target_dims, random_state=42, \n                            max_iter=200, init='nndsvda')\n        \n        # 2. NMF with different initialization\n        methods['NMF (Random)'] = NMF(n_components=self.target_dims, random_state=42, \n                                     max_iter=200, init='random')\n        \n        # 3. Dictionary Learning\n        methods['Dictionary Learning'] = DictionaryLearning(n_components=self.target_dims, \n                                                          alpha=0.1, random_state=42, \n                                                          max_iter=50)\n        \n        # Apply methods\n        for name, method in methods.items():\n            try:\n                print(f\"Applying {name}...\")\n                start_time = time.time()\n                \n                if 'NMF' in name:\n                    X_input = X_positive\n                else:\n                    X_input = X\n                \n                X_reduced = method.fit_transform(X_input)\n                \n                self.timings[name] = time.time() - start_time\n                self.results[name] = X_reduced\n                self.fitted_reducers[name] = method\n                self.metrics[name] = self.evaluate_reduction(X, X_reduced, y, name)\n                \n                print(f\"  {name}: {X_reduced.shape[1]} features, time: {self.timings[name]:.2f}s, \"\n                      f\"avg_corr: {self.metrics[name].get('avg_correlation', 0):.4f}\")\n                \n            except Exception as e:\n                print(f\"  {name} failed: {e}\")\n    \n    def apply_supervised_methods(self, X, y, feature_names):\n        \"\"\"Apply supervised dimensionality reduction methods\"\"\"\n        \n        print(\"\\n=== SUPERVISED METHODS ===\")\n        \n        # 1. PLS Regression\n        try:\n            print(\"Applying PLS...\")\n            start_time = time.time()\n            \n            pls = PLSRegression(n_components=min(self.target_dims, X.shape[1]//2))\n            X_reduced = pls.fit_transform(X, y)[0]  # Returns tuple, get X\n            \n            self.timings['PLS'] = time.time() - start_time\n            self.results['PLS'] = X_reduced\n            self.fitted_reducers['PLS'] = pls\n            self.metrics['PLS'] = self.evaluate_reduction(X, X_reduced, y, 'PLS')\n            \n            print(f\"  PLS: {X_reduced.shape[1]} features, time: {self.timings['PLS']:.2f}s, \"\n                  f\"avg_corr: {self.metrics['PLS'].get('avg_correlation', 0):.4f}\")\n            \n        except Exception as e:\n            print(f\"  PLS failed: {e}\")\n        \n        # 2. Linear Discriminant Analysis (for binned target)\n        try:\n            print(\"Applying LDA...\")\n            start_time = time.time()\n            \n            # Bin the target variable\n            y_binned = pd.qcut(y, q=10, labels=False, duplicates='drop')\n            n_classes = len(np.unique(y_binned))\n            \n            lda = LinearDiscriminantAnalysis(n_components=min(self.target_dims, n_classes-1))\n            X_reduced = lda.fit_transform(X, y_binned)\n            \n            self.timings['LDA'] = time.time() - start_time\n            self.results['LDA'] = X_reduced\n            self.fitted_reducers['LDA'] = lda\n            self.metrics['LDA'] = self.evaluate_reduction(X, X_reduced, y, 'LDA')\n            \n            print(f\"  LDA: {X_reduced.shape[1]} features, time: {self.timings['LDA']:.2f}s, \"\n                  f\"avg_corr: {self.metrics['LDA'].get('avg_correlation', 0):.4f}\")\n            \n        except Exception as e:\n            print(f\"  LDA failed: {e}\")\n    \n    def apply_ensemble_methods(self, X, y, feature_names):\n        \"\"\"Apply ensemble methods\"\"\"\n        \n        print(\"\\n=== ENSEMBLE METHODS ===\")\n        \n        if len(self.results) < 3:\n            print(\"Not enough methods for ensemble\")\n            return\n        \n        # 1. Two-stage reduction\n        print(\"Applying two-stage reduction...\")\n        try:\n            # Stage 1: Random Projection\n            rp = GaussianRandomProjection(n_components=min(200, X.shape[1]//2), random_state=42)\n            X_stage1 = rp.fit_transform(X)\n            \n            # Stage 2: PCA\n            pca_final = PCA(n_components=self.target_dims)\n            X_two_stage = pca_final.fit_transform(X_stage1)\n            \n            self.results['Two-Stage (RP+PCA)'] = X_two_stage\n            self.metrics['Two-Stage (RP+PCA)'] = self.evaluate_reduction(X, X_two_stage, y, 'Two-Stage (RP+PCA)')\n            \n            print(f\"  Two-Stage: {X_two_stage.shape[1]} features, \"\n                  f\"avg_corr: {self.metrics['Two-Stage (RP+PCA)'].get('avg_correlation', 0):.4f}\")\n            \n        except Exception as e:\n            print(f\"  Two-stage reduction failed: {e}\")\n        \n        # 2. Ensemble averaging - Fixed to handle different dimensionalities\n        print(\"Creating ensemble average...\")\n        try:\n            # Get top methods by average correlation\n            method_scores = [(m, self.metrics[m].get('avg_correlation', 0)) \n                           for m in self.results.keys() if self.metrics[m].get('avg_correlation', 0) > 0]\n            \n            if len(method_scores) >= 3:\n                # Group methods by their output dimensionality\n                dim_groups = {}\n                for method_name, score in method_scores:\n                    n_dims = self.results[method_name].shape[1]\n                    if n_dims not in dim_groups:\n                        dim_groups[n_dims] = []\n                    dim_groups[n_dims].append((method_name, score))\n                \n                # Find the most common dimensionality with at least 3 methods\n                valid_dims = [dim for dim, methods in dim_groups.items() if len(methods) >= 3]\n                \n                if valid_dims:\n                    # Choose the dimensionality closest to our target\n                    best_dim = min(valid_dims, key=lambda x: abs(x - self.target_dims))\n                    methods_same_dim = sorted(dim_groups[best_dim], key=lambda x: x[1], reverse=True)[:3]\n                    \n                    # Standardize and average only methods with same dimensionality\n                    embeddings = []\n                    for method_name, _ in methods_same_dim:\n                        X_method = StandardScaler().fit_transform(self.results[method_name])\n                        embeddings.append(X_method)\n                    \n                    # Now we can safely average since all have same shape\n                    X_ensemble = np.mean(embeddings, axis=0)\n                    \n                    self.results['Ensemble Average'] = X_ensemble\n                    self.metrics['Ensemble Average'] = self.evaluate_reduction(X, X_ensemble, y, 'Ensemble Average')\n                    \n                    print(f\"  Ensemble Average: {X_ensemble.shape[1]} features, \"\n                          f\"avg_corr: {self.metrics['Ensemble Average'].get('avg_correlation', 0):.4f}\")\n                    print(f\"    Based on: {[m[0] for m in methods_same_dim]}\")\n                else:\n                    print(\"  Not enough methods with same dimensionality for ensemble average\")\n            else:\n                print(\"  Not enough methods with positive correlation for ensemble\")\n            \n        except Exception as e:\n            print(f\"  Ensemble average failed: {e}\")\n    \n    def run_comprehensive_analysis(self, X, y, feature_names=None):\n        \"\"\"Run all dimensionality reduction methods\"\"\"\n        \n        if feature_names is None:\n            feature_names = [f'feature_{i}' for i in range(X.shape[1])]\n        \n        print(f\"\\nRunning comprehensive dimensionality reduction analysis\")\n        print(f\"Dataset shape: {X.shape}\")\n        print(f\"Target dimensions: {self.target_dims}\")\n        print(f\"Number of samples for analysis: {len(X)}\")\n        \n        # Apply all method categories\n        self.apply_linear_methods(X, y, feature_names)\n        self.apply_nonlinear_manifold_methods(X, y, feature_names)\n        self.apply_matrix_factorization_methods(X, y, feature_names)\n        self.apply_supervised_methods(X, y, feature_names)\n        self.apply_ensemble_methods(X, y, feature_names)\n        \n        return self.results, self.metrics\n    \n    def create_summary_report(self):\n        \"\"\"Create summary of all methods\"\"\"\n        \n        if not self.metrics:\n            print(\"No results to summarize\")\n            return None\n        \n        # Create summary dataframe\n        summary_data = []\n        \n        for method, metrics in self.metrics.items():\n            row = {\n                'Method': method,\n                'Avg_Correlation': metrics.get('avg_correlation', 0),\n                'Max_Correlation': metrics.get('max_correlation', 0),\n                'Std_Correlation': metrics.get('std_correlation', 0),\n                'N_Features': metrics.get('n_features', 0),\n                'Time_Seconds': self.timings.get(method, 0)\n            }\n            \n            # Add model metrics\n            for model in ['RandomForest', 'Ridge']:\n                row[f'{model}_Pearson'] = metrics.get(f'{model}_reduced_pearson', 0)\n            \n            summary_data.append(row)\n        \n        summary_df = pd.DataFrame(summary_data)\n        \n        # Remove methods with zero correlation\n        summary_df = summary_df[summary_df['Avg_Correlation'] > 0]\n        \n        if len(summary_df) == 0:\n            print(\"No successful methods to report\")\n            return None\n        \n        summary_df = summary_df.sort_values('Avg_Correlation', ascending=False)\n        \n        # Print summary\n        print(\"\\n\" + \"=\"*80)\n        print(\"DIMENSIONALITY REDUCTION SUMMARY\")\n        print(\"=\"*80)\n        \n        print(f\"\\nSuccessful methods: {len(summary_df)} out of {len(self.metrics)}\")\n        \n        print(\"\\nTop Methods by Average Correlation:\")\n        print(\"-\"*80)\n        print(f\"{'Method':<30} {'Avg Corr':>10} {'Max Corr':>10} {'Features':>10} {'Time (s)':>10}\")\n        print(\"-\"*80)\n        for idx, row in summary_df.head(15).iterrows():\n            print(f\"{row['Method']:<30} {row['Avg_Correlation']:>10.4f} \"\n                  f\"{row['Max_Correlation']:>10.4f} {int(row['N_Features']):>10} \"\n                  f\"{row['Time_Seconds']:>10.2f}\")\n        \n        # Model-specific performance\n        print(\"\\n\" + \"=\"*60)\n        print(\"MODEL-SPECIFIC PERFORMANCE\")\n        print(\"=\"*60)\n        \n        for model in ['RandomForest', 'Ridge']:\n            print(f\"\\nTop methods for {model}:\")\n            model_col = f'{model}_Pearson'\n            if model_col in summary_df.columns:\n                top_for_model = summary_df.nlargest(5, model_col)\n                for idx, row in top_for_model.iterrows():\n                    print(f\"  {row['Method']:<30} Pearson: {row[model_col]:.4f}\")\n        \n        return summary_df\n\n\ndef main():\n    \"\"\"Main execution function\"\"\"\n    \n    # Store start time\n    start_time = datetime.now()\n    \n    print(\"=\"*100)\n    print(\"COMPREHENSIVE DRW CRYPTO DIMENSIONALITY REDUCTION ANALYSIS\")\n    print(\"Enhanced Version with Larger Sample Sizes\")\n    print(\"=\"*100)\n    print(f\"Started at: {start_time}\")\n    \n    # Load data\n    try:\n        print(\"\\nLoading data...\")\n        # Try Kaggle path first\n        train_path = os.path.join(DATA_PATH, 'train.parquet')\n        df = pd.read_parquet(train_path)\n        \n        # Use larger subset for analysis\n        sample_size = min(DEFAULT_SAMPLE_SIZE, len(df))\n        if sample_size < len(df):\n            # Use last N samples (most recent data)\n            df_subset = df.tail(sample_size).copy()\n            print(f\"Loaded {len(df)} total samples, using last {sample_size} for analysis\")\n        else:\n            df_subset = df.copy()\n            print(f\"Using all {len(df)} samples for analysis\")\n        \n        # Separate features and target\n        feature_cols = [col for col in df_subset.columns if col not in ['timestamp', 'label']]\n        X = df_subset[feature_cols].values\n        y = df_subset['label'].values\n        \n        # Print data statistics\n        print(f\"\\nData statistics:\")\n        print(f\"  Features: {len(feature_cols)}\")\n        print(f\"  Target mean: {np.mean(y):.6f}\")\n        print(f\"  Target std: {np.std(y):.6f}\")\n        print(f\"  Target range: [{np.min(y):.6f}, {np.max(y):.6f}]\")\n        \n    except Exception as e:\n        print(f\"\\nError loading data: {e}\")\n        print(f\"Creating synthetic data for demonstration ({SYNTHETIC_SAMPLE_SIZE} samples)...\")\n        \n        np.random.seed(42)\n        n_samples = SYNTHETIC_SAMPLE_SIZE\n        \n        # Create feature names matching DRW dataset\n        feature_cols = [f'X{i}' for i in range(1, 891)]\n        feature_cols.extend(['bid_qty', 'ask_qty', 'buy_qty', 'sell_qty', 'volume'])\n        \n        # Generate synthetic data with some infinite values\n        X = np.random.randn(n_samples, len(feature_cols))\n        \n        # Add some infinite values to simulate real data\n        inf_mask = np.random.random(X.shape) < 0.001\n        X[inf_mask] = np.where(np.random.random(np.sum(inf_mask)) > 0.5, np.inf, -np.inf)\n        \n        # Make some features more important\n        for i, col in enumerate(feature_cols):\n            if col in IMPORTANT_FEATURES:\n                finite_mask = np.isfinite(X[:, i])\n                X[finite_mask, i] = X[finite_mask, i] * 2 + np.sin(np.arange(np.sum(finite_mask)) * 0.01)\n        \n        # Create target with correlation to important features\n        important_indices = [i for i, col in enumerate(feature_cols) if col in IMPORTANT_FEATURES][:5]\n        y = np.zeros(n_samples)\n        for idx in important_indices:\n            finite_mask = np.isfinite(X[:, idx])\n            y[finite_mask] += X[finite_mask, idx] * 0.1\n        y += np.random.randn(n_samples) * 0.05\n    \n    # Check for infinite values\n    n_inf = np.sum(np.isinf(X))\n    n_nan = np.sum(np.isnan(X))\n    print(f\"\\nOriginal data: {n_inf} infinite values, {n_nan} NaN values\")\n    print(f\"Data shape: {X.shape}\")\n    print(f\"Memory usage: {X.nbytes / 1024**2:.2f} MB\")\n    \n    # Enhanced preprocessing\n    print(\"\\n\" + \"=\"*80)\n    print(\"DATA PREPROCESSING\")\n    print(\"=\"*80)\n    \n    preprocessor = RobustDataPreprocessor()\n    X_processed, feature_names_filtered = preprocessor.fit_transform(X, y, feature_cols)\n    \n    print(f\"\\nPreprocessed data shape: {X_processed.shape}\")\n    print(f\"Features retained: {len(feature_names_filtered)}\")\n    print(f\"Memory usage after preprocessing: {X_processed.nbytes / 1024**2:.2f} MB\")\n    \n    # Run dimensionality reduction analysis\n    print(\"\\n\" + \"=\"*80)\n    print(\"DIMENSIONALITY REDUCTION ANALYSIS\")\n    print(\"=\"*80)\n    \n    # Test with multiple target dimensions\n    target_dimensions = [50, 100]\n    \n    all_results = {}\n    \n    for target_dim in target_dimensions:\n        print(f\"\\n{'='*60}\")\n        print(f\"ANALYSIS WITH TARGET DIMENSIONS = {target_dim}\")\n        print('='*60)\n        \n        reducer = ExtendedDimensionalityReducer(\n            target_dims=target_dim,\n            preserve_important_features=False  # Simplified for this demo\n        )\n        \n        results, metrics = reducer.run_comprehensive_analysis(\n            X_processed, y, feature_names_filtered\n        )\n        \n        # Create summary\n        summary_df = reducer.create_summary_report()\n        \n        all_results[target_dim] = {\n            'reducer': reducer,\n            'results': results,\n            'metrics': metrics,\n            'summary': summary_df\n        }\n    \n    # Visualization for best configuration\n    best_dim = 50\n    if best_dim in all_results and all_results[best_dim]['summary'] is not None:\n        summary_df = all_results[best_dim]['summary']\n        \n        if len(summary_df) > 0:\n            # Create visualizations\n            fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(16, 6))\n            \n            # Plot 1: Top methods bar chart\n            top_methods = summary_df.head(10)\n            ax1.barh(range(len(top_methods)), top_methods['Avg_Correlation'])\n            ax1.set_yticks(range(len(top_methods)))\n            ax1.set_yticklabels(top_methods['Method'])\n            ax1.set_xlabel('Average Correlation with Target')\n            ax1.set_title(f'Top 10 Dimensionality Reduction Methods (dims={best_dim})')\n            ax1.grid(True, alpha=0.3)\n            \n            # Plot 2: Time vs Performance scatter\n            ax2.scatter(summary_df['Time_Seconds'], summary_df['Avg_Correlation'], \n                       s=100, alpha=0.6)\n            \n            # Annotate top 5 methods\n            for idx, row in summary_df.head(5).iterrows():\n                ax2.annotate(row['Method'], \n                           (row['Time_Seconds'], row['Avg_Correlation']),\n                           xytext=(5, 5), textcoords='offset points', \n                           fontsize=8, alpha=0.8)\n            \n            ax2.set_xlabel('Computation Time (seconds)')\n            ax2.set_ylabel('Average Correlation with Target')\n            ax2.set_title('Speed vs Performance Trade-off')\n            ax2.grid(True, alpha=0.3)\n            \n            plt.tight_layout()\n            plt.show()\n    \n    # Final recommendations\n    print(\"\\n\" + \"=\"*80)\n    print(\"FINAL RECOMMENDATIONS FOR DRW CRYPTO COMPETITION\")\n    print(\"=\"*80)\n    \n    print(f\"\"\"\n    Based on analysis of {X_processed.shape[0]} samples:\n    \n    1. **DATA PREPROCESSING:**\n       - Winsorization effectively handles infinite values\n       - {len(feature_cols) - len(feature_names_filtered)} zero-variance features removed\n       - RobustScaler provides outlier-resistant scaling\n    \n    2. **TOP PERFORMING METHODS:**\"\"\")\n    \n    # Print top methods across all dimensions tested\n    for dim in target_dimensions:\n        if dim in all_results and all_results[dim]['summary'] is not None:\n            top_method = all_results[dim]['summary'].iloc[0]\n            print(f\"       - For {dim} dims: {top_method['Method']} \"\n                  f\"(corr: {top_method['Avg_Correlation']:.4f})\")\n    \n    print(\"\"\"\n    3. **COMPUTATIONAL CONSIDERATIONS:**\n       - Linear methods (PCA, Random Projection) scale well to large datasets\n       - Incremental PCA useful for very large datasets\n       - Two-stage reduction balances quality and speed\n    \n    4. **IMPLEMENTATION STRATEGY:**\n       - Use larger sample sizes for more robust estimates\n       - Consider ensemble methods for best performance\n       - Monitor memory usage with large datasets\n       - Save all preprocessing and reduction transformers\n    \"\"\")\n    \n    # Calculate total elapsed time\n    end_time = datetime.now()\n    elapsed_time = (end_time - start_time).total_seconds()\n    \n    print(f\"\\nTotal analysis time: {elapsed_time:.1f} seconds\")\n    print(f\"Completed at: {end_time}\")\n    print(\"\\n✅ Analysis complete with enhanced sample sizes!\")\n    \n    return all_results\n\n\nif __name__ == \"__main__\":\n    results = main()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-06-17T05:23:10.824553Z","iopub.execute_input":"2025-06-17T05:23:10.824873Z","execution_failed":"2025-06-17T12:40:09.999Z"}},"outputs":[],"execution_count":null}]}