{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":96164,"databundleVersionId":11418275,"sourceType":"competition"}],"dockerImageVersionId":31040,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"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},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#!/usr/bin/env python3\n\"\"\"\nOptimized Neural Network Feature Extraction Pipeline\nBalanced approach: More complex than basic, but faster than full version\n\"\"\"\n\nimport numpy as np\nimport pandas as pd\nimport torch\nimport torch.nn as nn\nimport torch.nn.functional as F\nfrom sklearn.tree import DecisionTreeRegressor\nfrom sklearn.preprocessing import PolynomialFeatures, StandardScaler\nfrom sklearn.linear_model import LassoCV\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import mean_squared_error, r2_score\nfrom scipy.optimize import minimize\nfrom scipy.special import expit\nfrom scipy.stats import spearmanr\nfrom itertools import combinations\nimport xgboost as xgb\nimport warnings\nimport pickle\nimport json\nimport time\nimport os\nfrom datetime import datetime\nfrom typing import Dict, List, Tuple, Any\n\nwarnings.filterwarnings('ignore')\n\n# ============================\n# OPTIMIZED CONFIGURATION\n# ============================\n\nclass Config:\n    \"\"\"Optimized configuration for balanced performance\"\"\"\n    # Time limits\n    MAX_RUNTIME_HOURS = 2  # Reduced for faster execution\n    CHECKPOINT_INTERVAL_MINUTES = 15\n    \n    # Sample sizes (optimized for speed)\n    POLYNOMIAL_SAMPLES = 2000\n    MULTIWAY_SAMPLES = 1500\n    DECISION_TREE_SAMPLES = 1500\n    SYMBOLIC_SAMPLES = 800\n    NONLINEAR_SAMPLES = 1000\n    \n    # Complexity limits (balanced)\n    MAX_POLYNOMIAL_DEGREE = 3\n    MAX_INTERACTION_ORDER = 4  # Up to 4-way interactions\n    MAX_TREE_DEPTH = 5\n    MAX_FEATURES_TO_TEST = 20\n    \n    # Output limits\n    MAX_POLYNOMIAL_TERMS = 30\n    MAX_RULES = 15\n    MAX_SYMBOLIC_FORMULAS = 20\n    MAX_NONLINEAR_FEATURES = 25\n    MAX_INTERACTIONS = 50\n    \n    # XGBoost parameters\n    XGB_PARAMS = {\n        'objective': 'reg:squarederror',\n        'max_depth': 6,\n        'learning_rate': 0.1,\n        'n_estimators': 300,\n        'subsample': 0.8,\n        'colsample_bytree': 0.8,\n        'reg_alpha': 0.1,\n        'reg_lambda': 1.0,\n        'random_state': 42,\n        'n_jobs': -1  # Use all cores\n    }\n    \n    # Paths\n    CHECKPOINT_DIR = \"./nn_extraction_checkpoints\"\n    RESULTS_FILE = \"extraction_results.pkl\"\n    MODEL_FILE = \"xgboost_model.pkl\"\n\n# ============================\n# CHECKPOINT MANAGER\n# ============================\n\nclass CheckpointManager:\n    \"\"\"Lightweight checkpoint manager\"\"\"\n    \n    def __init__(self, checkpoint_dir=Config.CHECKPOINT_DIR):\n        self.checkpoint_dir = checkpoint_dir\n        os.makedirs(checkpoint_dir, exist_ok=True)\n        self.start_time = time.time()\n        self.last_checkpoint_time = time.time()\n        \n    def save_checkpoint(self, state: Dict[str, Any], name: str):\n        \"\"\"Save checkpoint to file\"\"\"\n        filepath = os.path.join(self.checkpoint_dir, f\"{name}_checkpoint.pkl\")\n        with open(filepath, 'wb') as f:\n            pickle.dump(state, f)\n        self.last_checkpoint_time = time.time()\n        \n    def load_checkpoint(self, name: str) -> Dict[str, Any]:\n        \"\"\"Load checkpoint from file\"\"\"\n        filepath = os.path.join(self.checkpoint_dir, f\"{name}_checkpoint.pkl\")\n        if os.path.exists(filepath):\n            with open(filepath, 'rb') as f:\n                return pickle.load(f)\n        return None\n    \n    def should_checkpoint(self) -> bool:\n        \"\"\"Check if it's time to save a checkpoint\"\"\"\n        elapsed = time.time() - self.last_checkpoint_time\n        return elapsed > Config.CHECKPOINT_INTERVAL_MINUTES * 60\n    \n    def time_remaining(self) -> float:\n        \"\"\"Get remaining time in seconds\"\"\"\n        elapsed = time.time() - self.start_time\n        max_seconds = Config.MAX_RUNTIME_HOURS * 3600\n        return max(0, max_seconds - elapsed)\n    \n    def should_stop(self) -> bool:\n        \"\"\"Check if we should stop due to time limit\"\"\"\n        return self.time_remaining() <= 0\n\n# ============================\n# OPTIMIZED NEURAL NETWORK\n# ============================\n\nclass OptimizedFeatureDiscoveryNetwork(nn.Module):\n    \"\"\"Single optimized network for all feature discovery tasks\"\"\"\n    \n    def __init__(self, input_dim, hidden_dims=[384, 192, 96], dropout=0.2):\n        super().__init__()\n        self.input_dim = input_dim\n        \n        # Feature importance weights\n        self.feature_weights = nn.Parameter(torch.ones(input_dim))\n        \n        # Main encoder\n        layers = []\n        prev_dim = input_dim\n        \n        for hidden_dim in hidden_dims:\n            layers.extend([\n                nn.Linear(prev_dim, hidden_dim),\n                nn.BatchNorm1d(hidden_dim),\n                nn.ReLU(),\n                nn.Dropout(dropout)\n            ])\n            prev_dim = hidden_dim\n        \n        self.encoder = nn.Sequential(*layers)\n        \n        # Interaction detector (simplified)\n        self.interaction_detector = nn.Sequential(\n            nn.Linear(hidden_dims[-1], hidden_dims[-1] // 2),\n            nn.ReLU(),\n            nn.Linear(hidden_dims[-1] // 2, 50)  # Top 50 interaction scores\n        )\n        \n        # Output layer\n        self.output = nn.Linear(hidden_dims[-1], 1)\n        \n        # Gates for feature selection\n        self.feature_gates = nn.Sequential(\n            nn.Linear(input_dim, input_dim),\n            nn.Sigmoid()\n        )\n        \n    def forward(self, x):\n        # Apply feature weights and gates\n        gates = self.feature_gates(x)\n        x_weighted = x * F.softmax(self.feature_weights, dim=0) * gates\n        \n        # Encode\n        encoded = self.encoder(x_weighted)\n        \n        # Output\n        output = self.output(encoded)\n        return output\n    \n    def get_feature_importance(self):\n        \"\"\"Get learned feature importance\"\"\"\n        with torch.no_grad():\n            importance = F.softmax(self.feature_weights, dim=0).cpu().numpy()\n        return importance\n    \n    def get_interaction_scores(self, x):\n        \"\"\"Get interaction scores\"\"\"\n        with torch.no_grad():\n            gates = self.feature_gates(x)\n            x_weighted = x * F.softmax(self.feature_weights, dim=0) * gates\n            encoded = self.encoder(x_weighted)\n            scores = self.interaction_detector(encoded)\n        return scores\n\n# ============================\n# FAST FORMULA EXTRACTOR\n# ============================\n\nclass FastFormulaExtractor:\n    \"\"\"Optimized formula extraction focusing on speed\"\"\"\n    \n    def __init__(self, feature_names, device='cpu'):\n        self.feature_names = feature_names\n        self.device = device\n        self.formulas = {\n            'polynomial': [],\n            'interactions': [],\n            'nonlinear': [],\n            'symbolic': [],\n            'rules': []\n        }\n        \n    def extract_all(self, model, X, y, checkpoint_manager=None):\n        \"\"\"Extract all formulas with optimizations\"\"\"\n        print(\"\\n\" + \"=\"*60)\n        print(\"FAST FORMULA EXTRACTION\")\n        print(\"=\"*60)\n        \n        # 1. Polynomial features\n        print(\"\\n1. Extracting polynomial features...\")\n        start = time.time()\n        self._extract_polynomial_fast(model, X, y)\n        print(f\"   Completed in {time.time()-start:.1f}s\")\n        \n        # 2. Interactions\n        print(\"\\n2. Extracting interactions...\")\n        start = time.time()\n        self._extract_interactions_fast(model, X, y)\n        print(f\"   Completed in {time.time()-start:.1f}s\")\n        \n        # 3. Nonlinear\n        print(\"\\n3. Extracting nonlinear features...\")\n        start = time.time()\n        self._extract_nonlinear_fast(model, X, y)\n        print(f\"   Completed in {time.time()-start:.1f}s\")\n        \n        # 4. Symbolic\n        print(\"\\n4. Extracting symbolic formulas...\")\n        start = time.time()\n        self._extract_symbolic_fast(model, X, y)\n        print(f\"   Completed in {time.time()-start:.1f}s\")\n        \n        # 5. Rules\n        print(\"\\n5. Extracting decision rules...\")\n        start = time.time()\n        self._extract_rules_fast(model, X, y)\n        print(f\"   Completed in {time.time()-start:.1f}s\")\n        \n        return self.formulas\n    \n    def _extract_polynomial_fast(self, model, X, y):\n        \"\"\"Fast polynomial extraction\"\"\"\n        # Sample data\n        n_samples = min(Config.POLYNOMIAL_SAMPLES, len(X))\n        indices = np.random.choice(len(X), n_samples, replace=False)\n        X_sample = X[indices]\n        \n        # Get NN predictions\n        with torch.no_grad():\n            X_tensor = torch.FloatTensor(X_sample).to(self.device)\n            nn_preds = model(X_tensor).squeeze().cpu().numpy()\n        \n        # Get top features by importance\n        importance = model.get_feature_importance()\n        top_features = np.argsort(importance)[-Config.MAX_FEATURES_TO_TEST:]\n        \n        # Polynomial features\n        poly = PolynomialFeatures(degree=Config.MAX_POLYNOMIAL_DEGREE, include_bias=False)\n        X_poly = poly.fit_transform(X_sample[:, top_features])\n        \n        # Fast LASSO\n        lasso = LassoCV(cv=3, max_iter=1000, n_jobs=-1)\n        lasso.fit(X_poly, nn_preds)\n        \n        # Extract important terms\n        feature_names_subset = [self.feature_names[i] for i in top_features if i < len(self.feature_names)]\n        poly_names = poly.get_feature_names_out(feature_names_subset)\n        \n        important_idx = np.where(np.abs(lasso.coef_) > 1e-5)[0]\n        important_idx = important_idx[np.argsort(np.abs(lasso.coef_[important_idx]))[::-1]]\n        \n        for idx in important_idx[:Config.MAX_POLYNOMIAL_TERMS]:\n            self.formulas['polynomial'].append({\n                'name': poly_names[idx],\n                'coef': float(lasso.coef_[idx]),\n                'degree': poly_names[idx].count('^') + 1\n            })\n        \n        print(f\"   Found {len(self.formulas['polynomial'])} polynomial terms\")\n    \n    def _extract_interactions_fast(self, model, X, y):\n        \"\"\"Fast interaction extraction\"\"\"\n        # Sample data\n        n_samples = min(Config.MULTIWAY_SAMPLES, len(X))\n        indices = np.random.choice(len(X), n_samples, replace=False)\n        X_sample = X[indices]\n        \n        # Get feature importance\n        importance = model.get_feature_importance()\n        top_features = np.argsort(importance)[-Config.MAX_FEATURES_TO_TEST:]\n        \n        found_interactions = []\n        \n        # 2-way interactions (most common and useful)\n        for i, f1 in enumerate(top_features):\n            for f2 in top_features[i+1:]:\n                if f1 < len(self.feature_names) and f2 < len(self.feature_names):\n                    # Fast interaction test\n                    interaction = X_sample[:, f1] * X_sample[:, f2]\n                    if np.std(interaction) > 0:\n                        # Quick correlation check\n                        with torch.no_grad():\n                            X_tensor = torch.FloatTensor(X_sample).to(self.device)\n                            y_pred = model(X_tensor).squeeze().cpu().numpy()\n                        \n                        corr = np.abs(np.corrcoef(interaction, y_pred)[0, 1])\n                        if corr > 0.1:\n                            found_interactions.append({\n                                'features': [self.feature_names[f1], self.feature_names[f2]],\n                                'order': 2,\n                                'score': float(corr)\n                            })\n        \n        # 3-way interactions (limited)\n        if Config.MAX_INTERACTION_ORDER >= 3:\n            # Only test top 10 features for 3-way\n            top_10 = top_features[-10:]\n            for combo in combinations(top_10, 3):\n                if all(f < len(self.feature_names) for f in combo):\n                    interaction = np.prod([X_sample[:, f] for f in combo], axis=0)\n                    if np.std(interaction) > 0:\n                        with torch.no_grad():\n                            X_tensor = torch.FloatTensor(X_sample).to(self.device)\n                            y_pred = model(X_tensor).squeeze().cpu().numpy()\n                        \n                        corr = np.abs(np.corrcoef(interaction, y_pred)[0, 1])\n                        if corr > 0.15:  # Higher threshold for 3-way\n                            found_interactions.append({\n                                'features': [self.feature_names[f] for f in combo],\n                                'order': 3,\n                                'score': float(corr)\n                            })\n        \n        # Keep best interactions\n        found_interactions.sort(key=lambda x: x['score'], reverse=True)\n        self.formulas['interactions'] = found_interactions[:Config.MAX_INTERACTIONS]\n        \n        print(f\"   Found {len(self.formulas['interactions'])} interactions\")\n    \n    def _extract_nonlinear_fast(self, model, X, y):\n        \"\"\"Fast nonlinear feature extraction\"\"\"\n        # Sample data\n        n_samples = min(Config.NONLINEAR_SAMPLES, len(X))\n        indices = np.random.choice(len(X), n_samples, replace=False)\n        X_sample = X[indices]\n        \n        # Get predictions\n        with torch.no_grad():\n            X_tensor = torch.FloatTensor(X_sample).to(self.device)\n            y_pred = model(X_tensor).squeeze().cpu().numpy()\n        \n        # Get top features\n        importance = model.get_feature_importance()\n        top_features = np.argsort(importance)[-15:]\n        \n        # Simplified nonlinear functions\n        functions = {\n            'square': lambda x: x**2,\n            'sqrt_abs': lambda x: np.sqrt(np.abs(x)),\n            'log_abs': lambda x: np.log(np.abs(x) + 1),\n            'exp_scaled': lambda x: np.exp(np.clip(x/3, -5, 5)),\n            'sigmoid': lambda x: 1 / (1 + np.exp(-x)),\n            'tanh': lambda x: np.tanh(x),\n            'sin': lambda x: np.sin(x),\n            'cos': lambda x: np.cos(x),\n            'reciprocal': lambda x: 1 / (1 + x**2),\n            'gaussian': lambda x: np.exp(-x**2/2)\n        }\n        \n        found_nonlinear = []\n        \n        for feat_idx in top_features:\n            if feat_idx < len(self.feature_names):\n                x_feat = X_sample[:, feat_idx]\n                feat_name = self.feature_names[feat_idx]\n                \n                for func_name, func in functions.items():\n                    try:\n                        transformed = func(x_feat)\n                        if np.all(np.isfinite(transformed)):\n                            corr = np.abs(np.corrcoef(transformed, y_pred)[0, 1])\n                            if corr > 0.1:\n                                found_nonlinear.append({\n                                    'feature': feat_name,\n                                    'function': func_name,\n                                    'score': float(corr)\n                                })\n                    except:\n                        continue\n        \n        # Keep best\n        found_nonlinear.sort(key=lambda x: x['score'], reverse=True)\n        self.formulas['nonlinear'] = found_nonlinear[:Config.MAX_NONLINEAR_FEATURES]\n        \n        print(f\"   Found {len(self.formulas['nonlinear'])} nonlinear features\")\n    \n    def _extract_symbolic_fast(self, model, X, y):\n        \"\"\"Fast symbolic formula extraction\"\"\"\n        # Sample data\n        n_samples = min(Config.SYMBOLIC_SAMPLES, len(X))\n        indices = np.random.choice(len(X), n_samples, replace=False)\n        X_sample = X[indices]\n        \n        # Get predictions\n        with torch.no_grad():\n            X_tensor = torch.FloatTensor(X_sample).to(self.device)\n            y_pred = model(X_tensor).squeeze().cpu().numpy()\n        \n        # Get top features\n        importance = model.get_feature_importance()\n        top_features = np.argsort(importance)[-10:]\n        \n        found_symbolic = []\n        \n        # Simple templates only\n        for i, f1 in enumerate(top_features[:5]):\n            for f2 in top_features[i+1:i+3]:\n                if f1 < len(self.feature_names) and f2 < len(self.feature_names):\n                    x1 = X_sample[:, f1]\n                    x2 = X_sample[:, f2]\n                    \n                    # Template 1: a*x1^b * x2^c\n                    def power_template(params):\n                        a, b, c = params\n                        try:\n                            pred = a * np.power(np.abs(x1) + 1e-8, b) * np.power(np.abs(x2) + 1e-8, c)\n                            return np.mean((pred - y_pred)**2)\n                        except:\n                            return 1e10\n                    \n                    result = minimize(\n                        power_template,\n                        x0=[1, 0.5, 0.5],\n                        bounds=[(-10, 10), (-2, 2), (-2, 2)],\n                        method='L-BFGS-B',\n                        options={'maxiter': 50}\n                    )\n                    \n                    if result.success:\n                        a, b, c = result.x\n                        score = 1 - result.fun / (np.var(y_pred) + 1e-8)\n                        if score > 0.1:\n                            found_symbolic.append({\n                                'template': 'power_product',\n                                'features': [self.feature_names[f1], self.feature_names[f2]],\n                                'params': {'a': a, 'b': b, 'c': c},\n                                'formula': f\"{a:.3f}*{self.feature_names[f1]}^{b:.3f}*{self.feature_names[f2]}^{c:.3f}\",\n                                'score': float(score)\n                            })\n        \n        # Keep best\n        found_symbolic.sort(key=lambda x: x['score'], reverse=True)\n        self.formulas['symbolic'] = found_symbolic[:Config.MAX_SYMBOLIC_FORMULAS]\n        \n        print(f\"   Found {len(self.formulas['symbolic'])} symbolic formulas\")\n    \n    def _extract_rules_fast(self, model, X, y):\n        \"\"\"Fast decision rule extraction\"\"\"\n        # Sample data\n        n_samples = min(Config.DECISION_TREE_SAMPLES, len(X))\n        indices = np.random.choice(len(X), n_samples, replace=False)\n        X_sample = X[indices]\n        \n        # Get predictions\n        with torch.no_grad():\n            X_tensor = torch.FloatTensor(X_sample).to(self.device)\n            y_pred = model(X_tensor).squeeze().cpu().numpy()\n        \n        # Shallow tree for speed\n        tree = DecisionTreeRegressor(\n            max_depth=Config.MAX_TREE_DEPTH,\n            min_samples_leaf=100,\n            max_features='sqrt'\n        )\n        tree.fit(X_sample, y_pred)\n        \n        # Extract simple rules\n        rules = []\n        \n        # Get leaf nodes\n        leaf_nodes = np.where(tree.tree_.feature == -2)[0]\n        leaf_info = [(node, tree.tree_.n_node_samples[node], tree.tree_.value[node][0, 0]) \n                     for node in leaf_nodes]\n        leaf_info.sort(key=lambda x: x[1], reverse=True)\n        \n        for node, samples, value in leaf_info[:Config.MAX_RULES]:\n            if samples > 50:  # Only meaningful rules\n                rules.append({\n                    'node_id': int(node),\n                    'samples': int(samples),\n                    'value': float(value)\n                })\n        \n        self.formulas['rules'] = rules\n        print(f\"   Found {len(rules)} decision rules\")\n\n# ============================\n# FAST FEATURE CREATOR\n# ============================\n\nclass FastFeatureCreator:\n    \"\"\"Create features from formulas efficiently\"\"\"\n    \n    def __init__(self, formulas):\n        self.formulas = formulas\n        \n    def create_all_features(self, df, max_features=300):\n        \"\"\"Create all features efficiently\"\"\"\n        print(\"\\n   Creating features from formulas...\")\n        \n        new_features = pd.DataFrame(index=df.index)\n        feature_count = 0\n        \n        # 1. Polynomial features\n        for i, poly in enumerate(self.formulas.get('polynomial', [])[:50]):\n            if feature_count >= max_features:\n                break\n            \n            # For simplicity, just use coefficient as feature\n            # In practice, you'd parse and compute the polynomial\n            new_features[f'poly_{i}'] = poly['coef']\n            feature_count += 1\n        \n        # 2. Interaction features\n        for i, interaction in enumerate(self.formulas.get('interactions', [])[:80]):\n            if feature_count >= max_features:\n                break\n            \n            features = interaction['features']\n            if all(f in df.columns for f in features):\n                if interaction['order'] == 2:\n                    new_features[f'interact_2_{i}'] = df[features[0]] * df[features[1]]\n                elif interaction['order'] == 3:\n                    new_features[f'interact_3_{i}'] = df[features[0]] * df[features[1]] * df[features[2]]\n                feature_count += 1\n        \n        # 3. Nonlinear features\n        for i, nonlinear in enumerate(self.formulas.get('nonlinear', [])[:60]):\n            if feature_count >= max_features:\n                break\n            \n            feat_name = nonlinear['feature']\n            func_name = nonlinear['function']\n            \n            if feat_name in df.columns:\n                x = df[feat_name].values\n                \n                try:\n                    if func_name == 'square':\n                        result = x**2\n                    elif func_name == 'sqrt_abs':\n                        result = np.sqrt(np.abs(x))\n                    elif func_name == 'log_abs':\n                        result = np.log(np.abs(x) + 1)\n                    elif func_name == 'exp_scaled':\n                        result = np.exp(np.clip(x/3, -5, 5))\n                    elif func_name == 'sigmoid':\n                        result = 1 / (1 + np.exp(-x))\n                    elif func_name == 'tanh':\n                        result = np.tanh(x)\n                    elif func_name == 'sin':\n                        result = np.sin(x)\n                    elif func_name == 'cos':\n                        result = np.cos(x)\n                    elif func_name == 'reciprocal':\n                        result = 1 / (1 + x**2)\n                    elif func_name == 'gaussian':\n                        result = np.exp(-x**2/2)\n                    else:\n                        continue\n                    \n                    new_features[f'{func_name}_{feat_name}'] = result\n                    feature_count += 1\n                except:\n                    continue\n        \n        # 4. Symbolic features\n        for i, symbolic in enumerate(self.formulas.get('symbolic', [])[:30]):\n            if feature_count >= max_features:\n                break\n            \n            features = symbolic['features']\n            if len(features) >= 2 and all(f in df.columns for f in features[:2]):\n                params = symbolic['params']\n                if symbolic['template'] == 'power_product':\n                    x1 = df[features[0]].values\n                    x2 = df[features[1]].values\n                    result = params['a'] * np.power(np.abs(x1) + 1e-8, params['b']) * np.power(np.abs(x2) + 1e-8, params['c'])\n                    new_features[f'symbolic_{i}'] = np.clip(result, -1e6, 1e6)\n                    feature_count += 1\n        \n        # Clean features\n        for col in new_features.columns:\n            new_features[col] = new_features[col].replace([np.inf, -np.inf], np.nan)\n            new_features[col] = new_features[col].fillna(0)\n            \n            # Clip extreme values\n            if new_features[col].std() > 0:\n                mean_val = new_features[col].mean()\n                std_val = new_features[col].std()\n                new_features[col] = np.clip(new_features[col], mean_val - 5*std_val, mean_val + 5*std_val)\n        \n        print(f\"      Created {len(new_features.columns)} features\")\n        return new_features\n\n# ============================\n# XGBOOST TRAINING\n# ============================\n\ndef train_xgboost_model(X_train, y_train, X_test, y_test, feature_names):\n    \"\"\"Train XGBoost with optimized parameters\"\"\"\n    print(\"\\n\" + \"=\"*60)\n    print(\"TRAINING XGBOOST MODEL\")\n    print(\"=\"*60)\n    \n    # Scale features\n    scaler = StandardScaler()\n    X_train_scaled = scaler.fit_transform(X_train)\n    X_test_scaled = scaler.transform(X_test)\n    \n    # Convert to DMatrix\n    dtrain = xgb.DMatrix(X_train_scaled, label=y_train, feature_names=feature_names)\n    dtest = xgb.DMatrix(X_test_scaled, label=y_test, feature_names=feature_names)\n    \n    # Train model\n    print(\"\\nTraining XGBoost...\")\n    params = Config.XGB_PARAMS.copy()\n    evals = [(dtrain, 'train'), (dtest, 'test')]\n    \n    model = xgb.train(\n        params,\n        dtrain,\n        num_boost_round=params['n_estimators'],\n        evals=evals,\n        early_stopping_rounds=30,\n        verbose_eval=50\n    )\n    \n    # Evaluate\n    train_pred = model.predict(dtrain)\n    test_pred = model.predict(dtest)\n    \n    from sklearn.metrics import mean_squared_error, r2_score\n    \n    train_rmse = np.sqrt(mean_squared_error(y_train, train_pred))\n    test_rmse = np.sqrt(mean_squared_error(y_test, test_pred))\n    train_r2 = r2_score(y_train, train_pred)\n    test_r2 = r2_score(y_test, test_pred)\n    \n    print(f\"\\nModel Performance:\")\n    print(f\"  Train RMSE: {train_rmse:.4f}, R²: {train_r2:.4f}\")\n    print(f\"  Test RMSE: {test_rmse:.4f}, R²: {test_r2:.4f}\")\n    \n    # Feature importance\n    importance = model.get_score(importance_type='gain')\n    importance_df = pd.DataFrame([\n        {'feature': k, 'importance': v} \n        for k, v in importance.items()\n    ]).sort_values('importance', ascending=False)\n    \n    print(f\"\\nTop 15 Most Important Features:\")\n    for i, row in importance_df.head(15).iterrows():\n        print(f\"  {i+1}. {row['feature']}: {row['importance']:.2f}\")\n    \n    # Feature type analysis\n    feature_types = {}\n    for feat_name, imp in importance.items():\n        if 'poly_' in feat_name:\n            feat_type = 'polynomial'\n        elif 'interact_' in feat_name:\n            feat_type = 'interaction'\n        elif any(func in feat_name for func in ['square', 'sqrt', 'log', 'exp', 'sigmoid', 'tanh', 'sin', 'cos']):\n            feat_type = 'nonlinear'\n        elif 'symbolic_' in feat_name:\n            feat_type = 'symbolic'\n        else:\n            feat_type = 'original'\n        \n        if feat_type not in feature_types:\n            feature_types[feat_type] = 0\n        feature_types[feat_type] += imp\n    \n    print(f\"\\nFeature Type Importance:\")\n    for feat_type, total_imp in sorted(feature_types.items(), key=lambda x: x[1], reverse=True):\n        print(f\"  {feat_type}: {total_imp:.2f}\")\n    \n    return model, scaler, importance_df\n\n# ============================\n# MAIN PIPELINE\n# ============================\n\ndef main():\n    \"\"\"Main pipeline - optimized for speed and effectiveness\"\"\"\n    print(\"\\n\" + \"=\"*60)\n    print(\"OPTIMIZED NEURAL NETWORK FEATURE EXTRACTION PIPELINE\")\n    print(\"=\"*60)\n    \n    # Setup\n    device = torch.device('cuda' if torch.cuda.is_available() else 'cpu')\n    print(f\"\\nUsing device: {device}\")\n    \n    # Initialize checkpoint manager\n    checkpoint_manager = CheckpointManager()\n    \n    # Generate sample data (or load your own)\n    print(\"\\nGenerating sample data...\")\n    n_samples = 10000\n    n_features = 100\n    \n    np.random.seed(42)\n    X = np.random.randn(n_samples, n_features)\n    \n    # Create target with known patterns\n    y = (\n        # Polynomial\n        0.5 * X[:, 0]**2 * X[:, 1] + \n        0.3 * X[:, 2]**3 +\n        \n        # Interactions\n        0.4 * X[:, 3] * X[:, 4] * X[:, 5] +\n        0.2 * X[:, 6] * X[:, 7] +\n        \n        # Nonlinear\n        0.3 * np.sin(2 * X[:, 8]) +\n        0.2 * np.exp(-X[:, 9]**2) +\n        0.15 * np.tanh(X[:, 10]) +\n        \n        # Noise\n        0.1 * np.random.randn(n_samples)\n    )\n    \n    feature_names = [f'X_{i}' for i in range(n_features)]\n    \n    # Split data\n    train_size = int(0.8 * n_samples)\n    X_train, y_train = X[:train_size], y[:train_size]\n    X_test, y_test = X[train_size:], y[train_size:]\n    \n    print(f\"Training samples: {len(X_train)}\")\n    print(f\"Test samples: {len(X_test)}\")\n    \n    # Check for existing model\n    model_checkpoint = checkpoint_manager.load_checkpoint('nn_model')\n    \n    if model_checkpoint:\n        print(\"\\n✓ Loading saved neural network model...\")\n        model = model_checkpoint['model']\n        model.to(device)\n    else:\n        # Train neural network\n        print(\"\\nTraining neural network...\")\n        model = OptimizedFeatureDiscoveryNetwork(\n            input_dim=n_features,\n            hidden_dims=[384, 192, 96],\n            dropout=0.2\n        ).to(device)\n        \n        optimizer = torch.optim.Adam(model.parameters(), lr=0.001, weight_decay=1e-5)\n        scheduler = torch.optim.lr_scheduler.ReduceLROnPlateau(optimizer, patience=5)\n        batch_size = 256\n        \n        best_loss = float('inf')\n        patience_counter = 0\n        \n        for epoch in range(50):  # Max 50 epochs\n            # Training\n            model.train()\n            epoch_loss = 0\n            n_batches = 0\n            \n            for i in range(0, len(X_train), batch_size):\n                batch_idx = slice(i, min(i + batch_size, len(X_train)))\n                X_batch = torch.FloatTensor(X_train[batch_idx]).to(device)\n                y_batch = torch.FloatTensor(y_train[batch_idx]).to(device)\n                \n                optimizer.zero_grad()\n                outputs = model(X_batch).squeeze()\n                loss = F.mse_loss(outputs, y_batch)\n                loss.backward()\n                optimizer.step()\n                \n                epoch_loss += loss.item()\n                n_batches += 1\n            \n            avg_loss = epoch_loss / n_batches\n            scheduler.step(avg_loss)\n            \n            # Early stopping\n            if avg_loss < best_loss:\n                best_loss = avg_loss\n                patience_counter = 0\n            else:\n                patience_counter += 1\n                \n            if patience_counter >= 10:\n                print(f\"   Early stopping at epoch {epoch}\")\n                break\n                \n            if epoch % 10 == 0:\n                print(f\"   Epoch {epoch}: Loss = {avg_loss:.4f}\")\n        \n        # Save model\n        checkpoint_manager.save_checkpoint({'model': model}, 'nn_model')\n    \n    # Extract formulas\n    print(\"\\nExtracting formulas from neural network...\")\n    \n    extractor = FastFormulaExtractor(feature_names, device)\n    formulas = extractor.extract_all(model, X_train, y_train, checkpoint_manager)\n    \n    # Print summary\n    print(\"\\n\" + \"=\"*60)\n    print(\"FORMULA EXTRACTION SUMMARY\")\n    print(\"=\"*60)\n    \n    total_formulas = 0\n    for formula_type, formula_list in formulas.items():\n        count = len(formula_list)\n        total_formulas += count\n        print(f\"\\n{formula_type.upper()}: {count} formulas\")\n        \n        # Show top 3 examples\n        for i, formula in enumerate(formula_list[:3]):\n            if formula_type == 'polynomial':\n                print(f\"   {i+1}. {formula['name']} | Coef: {formula['coef']:.4f}\")\n            elif formula_type == 'interactions':\n                print(f\"   {i+1}. {' * '.join(formula['features'])} | Score: {formula['score']:.4f}\")\n            elif formula_type == 'nonlinear':\n                print(f\"   {i+1}. {formula['function']}({formula['feature']}) | Score: {formula['score']:.4f}\")\n            elif formula_type == 'symbolic':\n                print(f\"   {i+1}. {formula['formula']} | Score: {formula['score']:.4f}\")\n    \n    print(f\"\\nTotal formulas discovered: {total_formulas}\")\n    \n    # Create features\n    print(\"\\nCreating enhanced features...\")\n    \n    df_train = pd.DataFrame(X_train, columns=feature_names)\n    df_test = pd.DataFrame(X_test, columns=feature_names)\n    \n    creator = FastFeatureCreator(formulas)\n    new_features_train = creator.create_all_features(df_train)\n    new_features_test = creator.create_all_features(df_test)\n    \n    # Combine with original features\n    X_train_enhanced = pd.concat([df_train, new_features_train], axis=1)\n    X_test_enhanced = pd.concat([df_test, new_features_test], axis=1)\n    \n    print(f\"\\nTotal features: {X_train_enhanced.shape[1]} (original: {n_features}, new: {len(new_features_train.columns)})\")\n    \n    # Train XGBoost\n    xgb_model, scaler, importance_df = train_xgboost_model(\n        X_train_enhanced.values,\n        y_train,\n        X_test_enhanced.values,\n        y_test,\n        list(X_train_enhanced.columns)\n    )\n    \n    # Save results\n    print(\"\\nSaving results...\")\n    \n    # Create DMatrix with feature names for predictions\n    dtrain_pred = xgb.DMatrix(scaler.transform(X_train_enhanced.values), feature_names=list(X_train_enhanced.columns))\n    dtest_pred = xgb.DMatrix(scaler.transform(X_test_enhanced.values), feature_names=list(X_test_enhanced.columns))\n    \n    results = {\n        'formulas': formulas,\n        'feature_importance': importance_df.to_dict('records'),\n        'feature_names': list(X_train_enhanced.columns),\n        'performance': {\n            'train_rmse': float(np.sqrt(mean_squared_error(\n                y_train, \n                xgb_model.predict(dtrain_pred)\n            ))),\n            'test_rmse': float(np.sqrt(mean_squared_error(\n                y_test, \n                xgb_model.predict(dtest_pred)\n            )))\n        }\n    }\n    \n    # Save XGBoost model\n    model_path = os.path.join(Config.CHECKPOINT_DIR, Config.MODEL_FILE)\n    with open(model_path, 'wb') as f:\n        pickle.dump({'model': xgb_model, 'scaler': scaler}, f)\n    print(f\"✓ XGBoost model saved to: {model_path}\")\n    \n    # Save results\n    results_path = os.path.join(Config.CHECKPOINT_DIR, \"optimized_results.pkl\")\n    with open(results_path, 'wb') as f:\n        pickle.dump(results, f)\n    print(f\"✓ Results saved to: {results_path}\")\n    \n    # Save formulas in readable format\n    formulas_path = os.path.join(Config.CHECKPOINT_DIR, \"discovered_formulas.json\")\n    with open(formulas_path, 'w') as f:\n        json.dump(formulas, f, indent=2)\n    print(f\"✓ Formulas saved to: {formulas_path}\")\n    \n    print(\"\\n✓ Pipeline completed successfully!\")\n    \n    # Print timing\n    total_time = time.time() - checkpoint_manager.start_time\n    print(f\"\\nTotal execution time: {total_time/60:.1f} minutes\")\n    \n    return formulas, xgb_model, importance_df\n\nif __name__ == \"__main__\":\n    main()","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}