{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.7.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"gpu","dataSources":[{"sourceId":35332,"databundleVersionId":3723648,"sourceType":"competition"},{"sourceId":3739819,"sourceType":"datasetVersion","datasetId":2231132}],"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# 📆 Statement dates - Is there any valuable information for feature engineering?\n______________________________\n*updated 2022-07-20 [@roma-upgini](https://www.kaggle.com/romaupgini)*  🗣 Share this notebook: [Shareable Link](https://www.kaggle.com/romaupgini/statement-dates-to-use-or-not-to-use)\n\n\n## Four hypothesis on statement dates for feature engineering:\n\n1️⃣ There is last statement date seasonality for default rate and we can increase prediction accuracy by adding time dependent features for the last statement date (sin, cos, number of day, week etc). We have 31 days of statements (from 03/01 till 03/31). So, we can only check weekly seasonality.  \n\n2️⃣ There is an influence on default rate from holidays before/after the last statement, which might change income or spend structure for the households. We have to specify country for the holidays, but we can guess that for AMEX with a relatively high accuracy.  \n\n3️⃣ There is an influence on default rate from local economic situation, like emloyment, price inflation rates, central bank rates etc. Same - we'll need a country for that. However, 31 days of March 2018 most probably won't be enough to catch correlation with macroeconomical situation on training phase. But we have test set for April 2019 (public LB part) and October 2019 (private LB part) with 13-19 month shift from train period - enough for macroeconomics influence. It's worth to try.  \n\n4️⃣ Changes in statement dates for 13 month observation period has additional information correlated with default. From what I know about credit cards - consumer might change statement date, for example after salary date change. Which, in turn, might be the signal for employment change.\n\n**Let's check them one by one.**\n______________________________\n\n## Packages and functions\n\n📚 In this notebook we'll use:\n* [Upgini](https://github.com/upgini/upgini#readme) - Free data & feature enrichment library for machine learning: automatically enriches your training dataset with only the accuracy improving features from public and community shared data sources. <a href=\"https://github.com/upgini/upgini\">\n    <img src=\"https://img.shields.io/badge/GitHub-100000?style=for-the-badge&logo=github&logoColor=white\"  align='center'>\n</a>\n* [CuDF from RAPIDS.ai](https://github.com/rapidsai/cudf) - DataFrame on GPU\n\n**Switch on Internet and GPU for this kernel!**","metadata":{"id":"2vTwQk68nSeH"}},{"cell_type":"code","source":"%pip install -Uq upgini\nimport pandas as pd, numpy as np\nimport sklearn\nimport matplotlib.pyplot as plt, gc, os\nimport seaborn as sns\nimport cupy, cudf\n\n# RANDOM SEED\nSEED = 42\n# FILL NAN VALUE\nNAN_VALUE = -127\n\ndef read_file2cudf(path = '', usecols = None):\n    # LOAD DATAFRAME\n    if usecols is not None: df = cudf.read_parquet(path, columns=usecols)\n    else: df = cudf.read_parquet(path)\n    # REDUCE DTYPE FOR CUSTOMER AND DATE\n    df['customer_ID'] = df['customer_ID'].str[-16:].str.hex_to_int().astype('int64')\n    df.S_2 = cudf.to_datetime( df.S_2 )\n    print('shape of data:', df.shape)\n    return df\n\n# CALCULATE SIZE OF EACH SEPARATE TEST PART\ndef get_rows(customers, test, NUM_PARTS = 4, verbose = ''):\n    chunk = len(customers)//NUM_PARTS\n    if verbose != '':\n        print(f'We will process {verbose} data as {NUM_PARTS} separate parts.')\n        print(f'There will be {chunk} customers in each part (except the last part).')\n        print('Below are number of rows in each part:')\n    rows = []\n\n    for k in range(NUM_PARTS):\n        if k==NUM_PARTS-1: cc = customers[k*chunk:]\n        else: cc = customers[k*chunk:(k+1)*chunk]\n        s = test.loc[test.customer_ID.isin(cc)].shape[0]\n        rows.append(s)\n    if verbose != '': print( rows )\n    return rows,chunk\n\ndef xgb_amex(y_pred, y_true):\n    return 'amex', amex_metric_np(y_pred,y_true.get_label())\ndef lgb_amex_metric(y_pred, y_true):\n    return 'amex', amex_metric_np(y_pred,y_true.get_label()), True\n\n# code by @https://www.kaggle.com/yunchonggan\n# https://www.kaggle.com/competitions/amex-default-prediction/discussion/328020\ndef amex_metric_np(preds: np.ndarray, target: np.ndarray) -> float:\n    n_pos = np.sum(target)\n    n_neg = target.shape[0] - n_pos\n\n    indices = np.argsort(preds)[::-1]\n    preds, target = preds[indices], target[indices]\n\n    weight = 20.0 - target * 19.0\n    cum_norm_weight = (weight * (1 / weight.sum())).cumsum()\n    four_pct_mask = cum_norm_weight <= 0.04\n    d = np.sum(target[four_pct_mask]) / n_pos\n\n    lorentz = (target * (1 / n_pos)).cumsum()\n    gini = ((lorentz - cum_norm_weight) * weight).sum()\n\n    gini_max = 10 * n_neg * (1 - 19 / (n_pos + 20 * n_neg))\n\n    g = gini / gini_max\n    return 0.5 * (g + d)\n\n# official metric\nimport pandas as pd\ndef amex_metric(y_true: pd.DataFrame, y_pred: pd.DataFrame) -> float:\n\n    def top_four_percent_captured(y_true: pd.DataFrame, y_pred: pd.DataFrame) -> float:\n        df = (pd.concat([y_true, y_pred], axis='columns')\n              .sort_values('prediction', ascending=False))\n        df['weight'] = df['target'].apply(lambda x: 20 if x==0 else 1)\n        four_pct_cutoff = int(0.04 * df['weight'].sum())\n        df['weight_cumsum'] = df['weight'].cumsum()\n        df_cutoff = df.loc[df['weight_cumsum'] <= four_pct_cutoff]\n        return (df_cutoff['target'] == 1).sum() / (df['target'] == 1).sum()\n        \n    def weighted_gini(y_true: pd.DataFrame, y_pred: pd.DataFrame) -> float:\n        df = (pd.concat([y_true, y_pred], axis='columns')\n              .sort_values('prediction', ascending=False))\n        df['weight'] = df['target'].apply(lambda x: 20 if x==0 else 1)\n        df['random'] = (df['weight'] / df['weight'].sum()).cumsum()\n        total_pos = (df['target'] * df['weight']).sum()\n        df['cum_pos_found'] = (df['target'] * df['weight']).cumsum()\n        df['lorentz'] = df['cum_pos_found'] / total_pos\n        df['gini'] = (df['lorentz'] - df['random']) * df['weight']\n        return df['gini'].sum()\n\n    def normalized_weighted_gini(y_true: pd.DataFrame, y_pred: pd.DataFrame) -> float:\n        y_true_pred = y_true.rename(columns={'target': 'prediction'})\n        return weighted_gini(y_true, y_pred) / weighted_gini(y_true, y_true_pred)\n\n    g = normalized_weighted_gini(y_true, y_pred)\n    d = top_four_percent_captured(y_true, y_pred)\n\n    return 0.5 * (g + d)","metadata":{"_kg_hide-input":true,"id":"fRp2xBeInSeM","outputId":"4ab54a36-e31d-4cec-f054-d2269736bb10","execution":{"iopub.status.busy":"2022-07-24T17:43:40.479982Z","iopub.execute_input":"2022-07-24T17:43:40.480688Z","iopub.status.idle":"2022-07-24T17:44:02.239807Z","shell.execute_reply.started":"2022-07-24T17:43:40.480533Z","shell.execute_reply":"2022-07-24T17:44:02.238599Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## Quick data exploration - Last statement date","metadata":{"id":"VSP5UNEVnSeO"}},{"cell_type":"code","source":"train = pd.read_parquet('../input/amex-data-integer-dtypes-parquet-format/train.parquet', columns=['S_2','customer_ID'])\ntrain['S_2'] = pd.to_datetime(train['S_2'])\ntrain = train.groupby('customer_ID')['S_2'].agg('max').reset_index() #last statement only\n\ndf_train_labels = pd.read_csv('../input/amex-default-prediction/train_labels.csv')\ntrain = train.merge(df_train_labels, on='customer_ID')\ndel df_train_labels\n_ = gc.collect()\nprint(\"Shape of data: \",train.shape)","metadata":{"_kg_hide-input":true,"id":"aCfEBLomnSeP","execution":{"iopub.status.busy":"2022-07-24T17:44:02.242103Z","iopub.execute_input":"2022-07-24T17:44:02.242484Z","iopub.status.idle":"2022-07-24T17:44:08.396821Z","shell.execute_reply.started":"2022-07-24T17:44:02.242455Z","shell.execute_reply":"2022-07-24T17:44:08.395199Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"First, let's check number of statements and default rates by day. Where **default rate** is a ratio of customers with default on specific date.","metadata":{"id":"ttMeUMtJnSeP"}},{"cell_type":"code","source":"plot_df=train.groupby(\"S_2\").target.count()\nnew_train = pd.DataFrame(train.groupby(\"S_2\").target.mean())\nprint (f\"Number of days in train set: {new_train['target'].count()}\")\nprint (f\"Standard deviation of Default ratio: {new_train['target'].std()}\")\n\nfig, ax = plt.subplots(2,1,figsize = (20,8))\nplot_df.plot(title = \"# Statements\", ax = ax[0])\nnew_train.plot(title = \"% of Defaults\", ax = ax[1])\nplt.show()\ndel plot_df, ax, fig","metadata":{"_kg_hide-input":true,"id":"SqqMRmnwnSeQ","execution":{"iopub.status.busy":"2022-07-24T17:44:08.399391Z","iopub.execute_input":"2022-07-24T17:44:08.399903Z","iopub.status.idle":"2022-07-24T17:44:09.045124Z","shell.execute_reply.started":"2022-07-24T17:44:08.399839Z","shell.execute_reply":"2022-07-24T17:44:09.043791Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Interesting, most likely there are\n* either a strong card sales seasonality (most probably statement date derived from card activation date) - ie very few card sales on Sundays, and a big spike on Saturdays\n* or strong preferences on statement dates by customers themselfs\n\n>Disclaimer - I'm not an AMEX customer, but for most of the banks you can choose statement date, so I'm extrapolating here\n\nNext, default ratio by day.  \nThere is a deviation between days, but it's hard to guess weither it is significant or not.\n**And we have two options here:**\n1. Let's assume it's a time series, where **y** is a Default Ratio. Then, if there is an influence like in Hypothesis #3 (macroeconomic influence), it must be some trend component in this TS. So we can check Stationarity of TS using [Augmented Dickey-Fuller test](https://en.wikipedia.org/wiki/Augmented_Dickey–Fuller_test) with a significance level of less than 5%.  \nThe intuition behind this test is that it determines how strongly a time series is defined by a trend. However, there is a case, when it's doesn't help us - Hypothesis #1. As Augmented Dickey-Fuller test  won't detect seasonal component (it's a stationary TS). **So we have to use something different.**\n\n2. Let's auto generate **A LOT** of features from Holiday calendars, Workweek calendars, Political calendars, Sport calendars, sin, cos for month/week, add economic indicators and financial market data by using any data enrichment library. Then do feature selection with feature permutation. In this case - only features which has a statistically significant influence on model accuracy will be picked up. **Here we'll be able to test ALL Hypothesis #1, #2, #3 AT ONCE.**\n\nSo let's do Option 2, as it's quicker, using [Upgini](https://github.com/upgini/upgini#readme) - Free automated data and feature enrichment library for machine learning applications.   \nIt will add automatically a lot of external information about dates, holidays, events, financial markets, consumer sentiments, weather etc. all for the specific country / location. Than automatically checks for relevance (ie influence on prediction accuracy improvement) and select only features which will improve it.  \nFull [list of data scources and features, such as weather features, calendar features, financial features, etc ](https://github.com/upgini/upgini#-connected-data-sources-and-coverage)","metadata":{"id":"S2k-_pxQnSeQ"}},{"cell_type":"markdown","source":"## Hypothesis 1️⃣, 2️⃣, 3️⃣ test with Upgini automated data and feature enrichment library\n\nTo initiate search with Upgini library, you need to define so called [*search keys*](https://github.com/upgini/upgini#-search-key-types-we-support-more-is-coming) - a set of columns to join external data sources and features. In this competition we can use the following keys:\n\n1. Column **date** should be used as **SearchKey.DATE**.;  \n2. **Country** as \"US\" and \"UK\" (ISO-3166 country code), as most of AMEX customers are from US, next major market is UK.\n    \nWith this set of search keys, our X dataset will be matched with [different date-specific features](https://github.com/upgini/upgini#-connected-data-sources-and-coverage), taking into account the country. Than relevant selection and ranking will be done.  \n  \nTo start the search, we need to initiate *scikit-learn* compartible `FeaturesEnricher` transformer with appropriate **search** parameters.    \nAfter that, we can call the **fit** or **fit_transform**  method of `features_enricher`.","metadata":{"id":"iQiKzD_7nSeR"}},{"cell_type":"code","source":"from upgini import FeaturesEnricher, SearchKey\nfrom upgini.dataset import Dataset\n\nenricher = FeaturesEnricher(\n    date_format=\"%Y-%m-%d\",\n    search_keys={\"S_2\": SearchKey.DATE},\n    api_key = \"UOE9bdK62stiNVnZjQF6SeHccupKnFo78A1JkufO9Rw\",\n    country_code = \"US\", # change that to UK for another run\n)\nDataset.MIN_ROWS_COUNT = 20 #small X dataset, removed internal checks","metadata":{"id":"X9OIUSzGnSeS","execution":{"iopub.status.busy":"2022-07-24T17:44:09.047837Z","iopub.execute_input":"2022-07-24T17:44:09.048335Z","iopub.status.idle":"2022-07-24T17:44:10.781683Z","shell.execute_reply.started":"2022-07-24T17:44:09.048304Z","shell.execute_reply":"2022-07-24T17:44:10.780259Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"For `FeaturesEnricher.fit()` method, just like in all scikit-learn transformers, we should pass **X_train** as the first argument and **y_train** as the second argument.   \n**y_train** is needed to select **only relevant features & datasets, which will improve accuracy**. And rank new external features according to their prediction contribution, calculated as a SHAP values.  ","metadata":{"id":"MrAywH7qnSeT"}},{"cell_type":"code","source":"enricher.fit(\n    new_train.drop(columns=\"target\").reset_index(),\n    new_train[\"target\"]\n)\ndel enricher, train, new_train\n_ = gc.collect()","metadata":{"id":"Ne412cW5nSeT","execution":{"iopub.status.busy":"2022-07-24T17:44:10.783759Z","iopub.execute_input":"2022-07-24T17:44:10.785027Z","iopub.status.idle":"2022-07-24T17:45:44.995240Z","shell.execute_reply.started":"2022-07-24T17:44:10.784984Z","shell.execute_reply":"2022-07-24T17:45:44.993916Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 🏁 Conclusion for Hypothesis #1, #2 and #3\n\n**No relevant external features on dates found - both for US, and for UK.**  \nSo we have **to reject** Hypothesis #1,#2 and #3 for this training dataset.  \n\nIMHO, there must be some influence, but to catch that, we need different training data structured in a following way:\n\n1. 3 months of statements for default customers, minimum\n2. gap window between months in 3 month set - at least 6 months  \n\nThan we'll have 15 months observation window in a following schema: **1 + 6 months gap + 1 + 6 months gap + 1**","metadata":{"id":"bHBHJCc6nSeT"}},{"cell_type":"markdown","source":"## Hypothesis 4️⃣ test, using optimized public notebook from [@cdeotte](https://www.kaggle.com/code/cdeotte/xgboost-starter-0-793) (🙏)\n\nBaseline Public score for this notebook was **0.793**, local CV **0.794**    \nI made following impovements\n* CV with Stratification and 4 Folds\n* Changed CV metric from OOF Score to Average Score\n* Removed features with changes in cals methodology between train and LBs: B_29 and S_9\n* Added \"after-pay\" features\n* Removed NaN replacement\n* Permutation feature selection\n* More feature eng. with additional statistics on last observations (1500+ features)\n* Changed hyperparams for XGB\n\nPublic score after these changes **0.796** with local CV **0.79639** which is more consistent **\"Local CV to LB match\"** than before\n\nNow, let's calculate features from statement dates distances and compare results on the Public LB after enrichment with this new features.   \nFor **quick** estimation we'll do 2 steps:\n1. Calculate Distance between statement dates. Impute first observation / only one statement case with mean() value.\n2. Calculate Stat. features for Statement dates distance column\n\n### Quick feature engineering on distance between statement dates:","metadata":{"id":"JMg4lcscnSeU"}},{"cell_type":"code","source":"def feature_engineer(df):\n    cat_features = [\"B_30\",\"B_38\",\"D_114\",\"D_116\",\"D_117\",\"D_120\",\"D_126\",\"D_63\",\"D_64\",\"D_66\",\"D_68\"]\n    \n    # Initial feature selection to speed up fitting, based on @ambros\n    # https://www.kaggle.com/code/ambrosm/amex-lightgbm-quickstart/notebook\n    features_avg = ['B_1', 'B_2', 'B_3', 'B_4', 'B_5', 'B_6', 'B_8', 'B_9', 'B_10', 'B_11', 'B_12', 'B_13', 'B_14', 'B_15', 'B_16', 'B_17', 'B_18', 'B_19', 'B_20', 'B_21', 'B_22', 'B_23', 'B_24', 'B_25', 'B_28', 'B_30', 'B_32', 'B_33', 'B_37', 'B_38', 'B_39', 'B_40', 'B_41', 'B_42', 'D_39', 'D_41', 'D_42', 'D_43', 'D_44', 'D_45', 'D_46', 'D_47', 'D_48', 'D_50', 'D_51', 'D_53', 'D_54', 'D_55', 'D_58', 'D_59', 'D_60', 'D_61', 'D_62', 'D_65', 'D_66', 'D_69', 'D_70', 'D_71', 'D_72', 'D_73', 'D_74', 'D_75', 'D_76', 'D_77', 'D_78', 'D_80', 'D_82', 'D_84', 'D_86', 'D_91', 'D_92', 'D_94', 'D_96', 'D_103', 'D_104', 'D_108', 'D_112', 'D_113', 'D_114', 'D_115', 'D_117', 'D_118', 'D_119', 'D_120', 'D_121', 'D_122', 'D_123', 'D_124', 'D_125', 'D_126', 'D_128', 'D_129', 'D_131', 'D_132', 'D_133', 'D_134', 'D_135', 'D_136', 'D_140', 'D_141', 'D_142', 'D_144', 'D_145', 'P_2', 'P_3', 'P_4', 'R_1', 'R_2', 'R_3', 'R_7', 'R_8', 'R_9', 'R_10', 'R_11', 'R_14', 'R_15', 'R_16', 'R_17', 'R_20', 'R_21', 'R_22', 'R_24', 'R_26', 'R_27', 'S_3', 'S_5', 'S_6', 'S_7', 'S_11', 'S_12', 'S_13', 'S_15', 'S_16', 'S_18', 'S_22', 'S_23', 'S_25', 'S_26']\n    features_min = ['B_2', 'B_4', 'B_5', 'B_9', 'B_13', 'B_14', 'B_15', 'B_16', 'B_17', 'B_19', 'B_20', 'B_28', 'B_33', 'B_36', 'B_42', 'D_39', 'D_41', 'D_42', 'D_45', 'D_46', 'D_48', 'D_50', 'D_51', 'D_53', 'D_55', 'D_56', 'D_58', 'D_59', 'D_60', 'D_62', 'D_70', 'D_71', 'D_74', 'D_75', 'D_78', 'D_83', 'D_102', 'D_112', 'D_113', 'D_115', 'D_118', 'D_119', 'D_121', 'D_122', 'D_128', 'D_132', 'D_140', 'D_141', 'D_144', 'D_145', 'P_2', 'P_3', 'R_1', 'R_27', 'S_3', 'S_5', 'S_7', 'S_11', 'S_12', 'S_23', 'S_25']\n    features_max = ['B_1', 'B_2', 'B_3', 'B_4', 'B_5', 'B_6', 'B_7', 'B_8', 'B_9', 'B_10', 'B_12', 'B_13', 'B_14', 'B_15', 'B_16', 'B_17', 'B_18', 'B_19', 'B_21', 'B_23', 'B_24', 'B_25', 'B_30', 'B_33', 'B_37', 'B_38', 'B_39', 'B_40', 'B_42', 'D_39', 'D_41', 'D_42', 'D_43', 'D_44', 'D_45', 'D_46', 'D_47', 'D_48', 'D_49', 'D_50', 'D_52', 'D_55', 'D_56', 'D_58', 'D_59', 'D_60', 'D_61', 'D_63', 'D_64', 'D_65', 'D_70', 'D_71', 'D_72', 'D_73', 'D_74', 'D_76', 'D_77', 'D_78', 'D_80', 'D_82', 'D_84', 'D_91', 'D_102', 'D_105', 'D_107', 'D_110', 'D_111', 'D_112', 'D_115', 'D_116', 'D_117', 'D_118', 'D_119', 'D_121', 'D_122', 'D_123', 'D_124', 'D_125', 'D_126', 'D_128', 'D_131', 'D_132', 'D_133', 'D_134', 'D_135', 'D_136', 'D_138', 'D_140', 'D_141', 'D_142', 'D_144', 'D_145', 'P_2', 'P_3', 'P_4', 'R_1', 'R_3', 'R_5', 'R_6', 'R_7', 'R_8', 'R_10', 'R_11', 'R_14', 'R_17', 'R_20', 'R_26', 'R_27', 'S_3', 'S_5', 'S_7', 'S_8', 'S_11', 'S_12', 'S_13', 'S_15', 'S_16', 'S_22', 'S_23', 'S_24', 'S_25', 'S_26', 'S_27']\n    features_last = ['B_1', 'B_2', 'B_3', 'B_4', 'B_5', 'B_6', 'B_7', 'B_8', 'B_9', 'B_10', 'B_11', 'B_12', 'B_13', 'B_14', 'B_15', 'B_16', 'B_17', 'B_18', 'B_19', 'B_20', 'B_21', 'B_22', 'B_23', 'B_24', 'B_25', 'B_26', 'B_28', 'B_30', 'B_32', 'B_33', 'B_36', 'B_37', 'B_38', 'B_39', 'B_40', 'B_41', 'B_42', 'D_39', 'D_41', 'D_42', 'D_43', 'D_44', 'D_45', 'D_46', 'D_47', 'D_48', 'D_49', 'D_50', 'D_51', 'D_52', 'D_53', 'D_54', 'D_55', 'D_56', 'D_58', 'D_59', 'D_60', 'D_61', 'D_62', 'D_63', 'D_64', 'D_65', 'D_69', 'D_70', 'D_71', 'D_72', 'D_73', 'D_75', 'D_76', 'D_77', 'D_78', 'D_79', 'D_80', 'D_81', 'D_82', 'D_83', 'D_86', 'D_91', 'D_96', 'D_105', 'D_106', 'D_112', 'D_114', 'D_119', 'D_120', 'D_121', 'D_122', 'D_124', 'D_125', 'D_126', 'D_127', 'D_130', 'D_131', 'D_132', 'D_133', 'D_134', 'D_138', 'D_140', 'D_141', 'D_142', 'D_145', 'P_2', 'P_3', 'P_4', 'R_1', 'R_2', 'R_3', 'R_4', 'R_5', 'R_6', 'R_7', 'R_8', 'R_9', 'R_10', 'R_11', 'R_12', 'R_13', 'R_14', 'R_15', 'R_19', 'R_20', 'R_26', 'R_27', 'S_3', 'S_5', 'S_6', 'S_7', 'S_8', 'S_11', 'S_12', 'S_13', 'S_16', 'S_19', 'S_20', 'S_22', 'S_23', 'S_24', 'S_25', 'S_26', 'S_27']\n    features_last = list(set(features_last)-set(cat_features))\n    features_max = list(set(features_max)-set(cat_features))\n    features_min = list(set(features_min)-set(cat_features))\n    features_avg = list(set(features_avg)-set(cat_features))\n    \n    # Drop non stable features for train-test, based on % of NaNs\n    #https://www.kaggle.com/code/onodera1/amex-eda-comparison-of-training-and-test-data    \n    df.drop([\"B_29\",\"S_9\"], axis=1, inplace = True)\n    \n    # Hypothesis #4 - retrieve info from statement dates as distance between the dates\n    # Than calculate 'mean', 'std', 'max', 'last' statistics for distances\n    # cudf doesn't support diff() as GroupBy function, slow pandas DF used\n    temp = df[[\"customer_ID\",\"S_2\"]].to_pandas()\n    temp[\"SDist\"]=temp.groupby(\"customer_ID\")[\"S_2\"].diff() / np.timedelta64(1, 'D')\n    # Impute with average distance 30.53 days\n    temp['SDist'].fillna(30.53, inplace=True)\n    df = cudf.concat([df,cudf.from_pandas(temp[\"SDist\"])], axis=1)\n    del temp\n    _ = gc.collect()\n    features_last.append('SDist')\n    features_avg.append('SDist')\n    features_max.append('SDist')\n    features_min.append('SDist')\n    \n    #https://www.kaggle.com/competitions/amex-default-prediction/discussion/328514\n    df.loc[(df.R_13==0) & (df.R_17==0) & (df.R_20==0) & (df.R_8==0), 'R_6'] = 0\n    df.loc[df.B_39==-1, 'B_36'] = 0\n    \n    # Compute \"after pay\" features\n    # https://www.kaggle.com/code/jiweiliu/rapids-cudf-feature-engineering-xgb\n    for bcol in [f'B_{i}' for i in [11,14,17]]+['D_39','D_131']+[f'S_{i}' for i in [16,23]]:\n        for pcol in ['P_2','P_3']:\n            if bcol in df.columns:\n                df[[f'{bcol}-{pcol}']] = df[bcol] - df[pcol]\n                features_last.append(f'{bcol}-{pcol}')\n                features_avg.append(f'{bcol}-{pcol}')\n                features_max.append(f'{bcol}-{pcol}')\n                features_min.append(f'{bcol}-{pcol}')\n                \n    # BASIC FEATURE ENGINEERING\n    # https://www.kaggle.com/code/huseyincot/amex-agg-data-how-it-created\n    # https://www.kaggle.com/code/jiweiliu/rapids-cudf-feature-engineering-xgb\n    \n    test_num_last = df.groupby(\"customer_ID\")[features_last].agg(['last','first'])\n    test_num_last.columns = ['_'.join(x) for x in test_num_last.columns]\n    test_num_min = df.groupby(\"customer_ID\")[features_min].agg(['min'])\n    test_num_min.columns = ['_'.join(x) for x in test_num_min.columns]\n    test_num_max = df.groupby(\"customer_ID\")[features_max].agg(['max'])\n    test_num_max.columns = ['_'.join(x) for x in test_num_max.columns]\n    test_num_avg = df.groupby(\"customer_ID\")[features_avg].agg(['mean'])\n    test_num_avg.columns = ['_'.join(x) for x in test_num_avg.columns]\n    test_num_std = df.groupby(\"customer_ID\")[list(set().union(features_avg,features_last,features_min,features_max))].agg(['std','quantile'])\n    test_num_std.columns = ['_'.join(x) for x in test_num_std.columns]\n\n    test_cat_agg = df.groupby(\"customer_ID\")[cat_features].agg(['last','first'])\n    test_cat_agg.columns = ['_'.join(x) for x in test_cat_agg.columns]\n   \n    #add last statement date, statements count and \"new customer\" category (LT=0.5)\n    test_date_agg = df.groupby(\"customer_ID\")[[\"S_2\",\"B_3\",\"D_104\"]].agg(['last','count'])\n    test_date_agg.columns = ['_'.join(x) for x in test_date_agg.columns]\n    test_date_agg.rename(columns = {'S_2_count':'LT','S_2_last':'S_2'}, inplace = True)\n    test_date_agg.loc[(test_date_agg.B_3_last.isnull()) & (test_date_agg.LT==1),'LT'] = 0.5\n    test_date_agg.loc[(test_date_agg.D_104_last.isnull()) & (test_date_agg.LT==1),'LT'] = 0.5\n    test_date_agg.drop([\"B_3_last\",\"D_104_last\",\"B_3_count\",\"D_104_count\"], axis=1, inplace = True)\n    \n    df = cudf.concat([test_date_agg, test_num_last, test_num_min, test_num_max, test_num_avg, test_num_std, test_cat_agg], axis=1)\n    del test_date_agg, test_num_last, test_num_min, test_num_max, test_num_avg, test_num_std, test_cat_agg\n    \n    # Ratios/diffs on last values as features, based on @ragnar123\n    # https://www.kaggle.com/code/ragnar123/amex-lgbm-dart-cv-0-7977\n    for col in list(set().union(features_last,features_avg)):\n        try:\n            df[f'{col}_last_first_div'] = df[f'{col}_last'] / df[f'{col}_first']\n            df[f'{col}_last_mean_sub'] = df[f'{col}_last'] - df[f'{col}_mean']\n            df[f'{col}_last_mean_div'] = df[f'{col}_last'] / df[f'{col}_mean']\n            df[f'{col}_last_max_div'] = df[f'{col}_last'] / df[f'{col}_max']\n            df[f'{col}_last_min_div'] = df[f'{col}_last'] / df[f'{col}_min']\n        except:\n            pass\n        \n    print('shape after engineering', df.shape )\n    return df","metadata":{"execution":{"iopub.status.busy":"2022-07-20T12:19:39.948950Z","iopub.execute_input":"2022-07-20T12:19:39.949774Z","iopub.status.idle":"2022-07-20T12:19:40.007956Z","shell.execute_reply.started":"2022-07-20T12:19:39.949736Z","shell.execute_reply":"2022-07-20T12:19:40.006654Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"train = []\n_ = gc.collect()\n# raddar Kaggle dataset\n# https://www.kaggle.com/datasets/raddar/amex-data-integer-dtypes-parquet-format\nPATH=\"../input/amex-data-integer-dtypes-parquet-format/train.parquet\"\ntrain = read_file2cudf(path = PATH)\ntrain = feature_engineer(train)\n\n# ADD TARGETS\ntargets = cudf.read_csv('../input/amex-default-prediction/train_labels.csv')\n#targets = cudf.read_csv('/content/drive/MyDrive/colab/train_labels.csv')\ntargets['customer_ID'] = targets['customer_ID'].str[-16:].str.hex_to_int().astype('int64')\ntargets = targets.set_index('customer_ID')\ntrain = train.merge(targets, left_index=True, right_index=True, how='left')\ntrain.target = train.target.astype('int8')\ndel targets\n\n# cudf merge above randomly shuffles rows\ntrain = train.sort_index().reset_index()\n\n# FEATURES\n# remove S_2 from FEATURES list\nFEATURES = train.columns[2:-1]\nprint(f'There are {len(FEATURES)} features!')","metadata":{"_kg_hide-input":true,"id":"ph7XZYHJnSeV","outputId":"17b72947-9f7f-4ff7-c6c9-63d300d92bda","execution":{"iopub.status.busy":"2022-07-20T12:19:40.027764Z","iopub.execute_input":"2022-07-20T12:19:40.028330Z","iopub.status.idle":"2022-07-20T12:22:31.075638Z","shell.execute_reply.started":"2022-07-20T12:19:40.028293Z","shell.execute_reply":"2022-07-20T12:22:31.074494Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### XGB Training on GPU\nTraining with 4 folds, it will take 1h 30min on Kaggle's P100 GPU","metadata":{"id":"QxnyOolrnSeV"}},{"cell_type":"code","source":"%%time\nfrom sklearn.model_selection import StratifiedKFold\nimport xgboost as xgb\n\n# FOLDS PER MODEL, as number of weeks in a month\nFOLDS = 4\n\n# XGB MODEL PARAMETERS\nxgb_parms = {\n            'objective': 'binary:logitraw', \n            'tree_method': 'gpu_hist',\n            'predictor':'gpu_predictor',\n            'max_depth': 7,\n            'subsample':0.88,\n            'colsample_bytree': 0.1,\n            'gamma':1.5,\n            'min_child_weight':8,\n            'lambda': 50,\n            'eta':0.03,\n            'learning_rate':0.02,\n            'random_state':SEED\n    }\n\nimportances = []\noof = []\nTRAIN_SUBSAMPLE = 1.0\n_ = gc.collect()\n\nfor j, cat_index in enumerate ([train[train.LT>0].index]):\n    score = 0\n    skf = StratifiedKFold(n_splits=FOLDS, shuffle = True, random_state=SEED)\n    for fold,(train_idx, valid_idx) in enumerate(skf.split(\n                train.loc[cat_index,:][[\"customer_ID\"]],\n                train.loc[cat_index,:].target.values.get())):\n\n        # TRAIN WITH SUBSAMPLE OF TRAIN FOLD DATA\n        if TRAIN_SUBSAMPLE<1.0:\n            np.random.seed(SEED)\n            train_idx = np.random.choice(train_idx, \n                           int(len(train_idx)*TRAIN_SUBSAMPLE), replace=False)\n            np.random.seed(None)\n\n        print('#'*25)\n        print('### Fold',fold+1)\n        print('### Train size',len(train_idx),'Valid size',len(valid_idx))\n        print(f'### Training with {int(TRAIN_SUBSAMPLE*100)}% fold data...')\n        print('#'*25)\n\n        # TRAIN, VALID, TEST FOR FOLD K\n        y_valid = train.loc[valid_idx, 'target']\n        dtrain = xgb.DMatrix(data=train.loc[train_idx, FEATURES],\n                             label=train.loc[train_idx, 'target'],\n                             )\n        dvalid = xgb.DMatrix(data=train.loc[valid_idx, FEATURES],\n                             label=y_valid,\n                             )\n        \n        # TRAIN MODEL FOLD K\n        model = xgb.train(xgb_parms, \n                    dtrain=dtrain,\n                    evals=[(dtrain,'train'),(dvalid,'valid')],\n                    num_boost_round=8000,\n                    early_stopping_rounds=1800,\n                    custom_metric=xgb_amex,\n                    maximize=True,\n                    verbose_eval=200) \n        model.save_model(f'XGB_fold{fold}_LT{j}.json')\n        del dtrain\n        _ = gc.collect()\n\n        # GET FEATURE IMPORTANCE FOR FOLD K\n        dd = model.get_score(importance_type='weight')\n        df = pd.DataFrame({'feature':dd.keys(),f'importance_{fold}':dd.values()})\n        importances.append(df)\n\n        # INFER OOF FOLD K\n        oof_preds = model.predict(dvalid, iteration_range=(0,model.best_ntree_limit))\n        acc = amex_metric(pd.DataFrame({'target':y_valid.values.get()}), \n                                        pd.DataFrame({'prediction':oof_preds}))\n        print('Kaggle Metric =',acc,'\\n')\n        score += acc\n\n        # SAVE OOF\n        df = train.loc[valid_idx, ['customer_ID','target'] ].to_pandas()\n        df['oof_pred'] = oof_preds\n        oof.append( df )\n\n        del dvalid, y_valid, model, dd, df\n        _ = gc.collect()\n\n    score /= FOLDS\n    print('Average CV Kaggle Metric for group =',score)\n\nprint('#'*25)\noof = pd.concat(oof,axis=0,ignore_index=True).set_index('customer_ID')\nscore = amex_metric(pd.DataFrame({'target':oof.target.values}), \n                                pd.DataFrame({'prediction':oof.oof_pred.values}))\nprint('OOF CV Kaggle Metric =',score)\n# CLEAN RAM\ndel oof, skf, cat_index\ndel train\n_ = gc.collect()","metadata":{"_kg_hide-input":true,"scrolled":true,"id":"Cg4IWfUfnSeV","outputId":"717192b8-2389-4497-8d29-7432e717ceb4","execution":{"iopub.status.busy":"2022-07-20T12:22:31.077256Z","iopub.execute_input":"2022-07-20T12:22:31.077893Z","iopub.status.idle":"2022-07-20T13:45:14.837208Z","shell.execute_reply.started":"2022-07-20T12:22:31.077856Z","shell.execute_reply":"2022-07-20T13:45:14.836320Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Local CV with the new features has a score **0.79601**  \nBaseline solution without Statement distance features had **0.79639** on local CV.   \nKeep going.\n\n### Feature importance of Statement Dates features\nLet's check feature importance for TOP 20 vars and for \"SDist\" vars (derived from Statement dates distance)","metadata":{"id":"aB4jyisCnSeW"}},{"cell_type":"code","source":"import matplotlib.pyplot as plt\n\ndf = importances[0].copy()\nfor k in range(1,FOLDS*(j+1)): df = df.merge(importances[k], on='feature', how='left')\ndf['importance'] = df.iloc[:,1:].mean(axis=1)\ndf = df.sort_values('importance',ascending=False)\n\nNUM_FEATURES = 40\nplt.figure(figsize=(10,5*NUM_FEATURES//10))\nplt.barh(np.arange(NUM_FEATURES,0,-1), df.importance.values[:NUM_FEATURES])\nplt.yticks(np.arange(NUM_FEATURES,0,-1), df.feature.values[:NUM_FEATURES])\nplt.title(f'Feature Importance - Top {NUM_FEATURES}')\nplt.show()\n\ndf = df[df.feature.str.find(\"SDis\") != -1]\nplt.figure(figsize=(10,5*df.shape[0]//10))\nplt.barh(np.arange(df.shape[0],0,-1), df.importance.values[:df.shape[0]])\nplt.yticks(np.arange(df.shape[0],0,-1), df.feature.values[:df.shape[0]])\nplt.title(f'XGB Feature Importance - Statement date features')\nplt.show()\n\ndel df, plt, importances\n_ = gc.collect()","metadata":{"_kg_hide-input":true,"id":"LSOGxNTvnSeW","outputId":"f9b07303-f373-423b-f7e5-26a2b5e5bf69","execution":{"iopub.status.busy":"2022-07-20T13:46:43.549492Z","iopub.execute_input":"2022-07-20T13:46:43.550541Z","iopub.status.idle":"2022-07-20T13:46:44.323570Z","shell.execute_reply.started":"2022-07-20T13:46:43.550476Z","shell.execute_reply":"2022-07-20T13:46:44.322559Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"**11 features** from Statement distance was selected, but none of them in TOP 40.   \nSo we might notice small improvement from them on Public LB.  \nLet's check that.  \n\n### Process Test Data, Predict and Submit\nWe will load @raddar dataset from [here][1] with discussion [here][2].\n\n[1]: https://www.kaggle.com/datasets/raddar/amex-data-integer-dtypes-parquet-format\n[2]: https://www.kaggle.com/competitions/amex-default-prediction/discussion/328514","metadata":{"id":"kBZZvKyBnSeW"}},{"cell_type":"code","source":"# COMPUTE SIZE OF 4 PARTS FOR TEST DATA\nNUM_PARTS = 4\nTEST_PATH = '../input/amex-data-integer-dtypes-parquet-format/test.parquet'\n\nprint(f'Reading test data...')\ntest = read_file2cudf(path = TEST_PATH, usecols = ['customer_ID','S_2'])\ncustomers = test[['customer_ID']].drop_duplicates().sort_index().values.flatten()\nrows,num_cust = get_rows(customers, test[['customer_ID']], NUM_PARTS = NUM_PARTS, verbose = 'test')\n\n# INFER TEST DATA IN PARTS\nskip_rows = 0\nskip_cust = 0\ntest_preds = []\n\nfor k in range(NUM_PARTS):\n    \n    # READ PART OF TEST DATA\n    print(f'\\nReading test data...')\n    test = read_file2cudf(path = TEST_PATH)\n    test = test.iloc[skip_rows:skip_rows+rows[k]]\n    skip_rows += rows[k]\n    print(f'=> Test part {k+1} has shape', test.shape)\n    \n    # PROCESS AND FEATURE ENGINEER PART OF TEST DATA\n    test = feature_engineer(test)\n    if k==NUM_PARTS-1: test = test.loc[customers[skip_cust:]]\n    else: test = test.loc[customers[skip_cust:skip_cust+num_cust]]\n    skip_cust += num_cust\n    \n    for j, cat_index in enumerate ([test[test.LT>0].index]):\n        # XGB\n        dtest = xgb.DMatrix(data=test.loc[cat_index,:][FEATURES])\n        model = xgb.Booster()\n        model.load_model(f'XGB_fold0_LT{j}.json')\n        preds = model.predict(dtest, iteration_range=(0,model.best_ntree_limit))\n        for f in range(1,FOLDS):\n            model.load_model(f'XGB_fold{f}_LT{j}.json')\n            preds += model.predict(dtest, iteration_range=(0,model.best_ntree_limit))\n        del dtest, model\n        _ = gc.collect()\n        preds /= FOLDS\n        # SAVE\n        df =  test.loc[cat_index].reset_index()[[\"customer_ID\"]].to_pandas()\n        df['prediction'] = preds\n        test_preds.append(df)\n        del df, preds\n        _ = gc.collect()\n\n    # CLEAN MEMORY\n    del test\n    _ = gc.collect()","metadata":{"_kg_hide-input":true,"id":"uXl-SddunSeW","outputId":"e240a589-2fec-4130-e194-4bbda97520fe","execution":{"iopub.status.busy":"2022-07-20T13:46:56.935663Z","iopub.execute_input":"2022-07-20T13:46:56.936659Z","iopub.status.idle":"2022-07-20T13:58:12.168324Z","shell.execute_reply.started":"2022-07-20T13:46:56.936607Z","shell.execute_reply":"2022-07-20T13:58:12.167336Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# WRITE SUBMISSION FILE\ntest = cudf.DataFrame.from_pandas(pd.concat(test_preds,axis=0,ignore_index=True).set_index(\"customer_ID\"))\nsub = cudf.read_csv('../input/amex-default-prediction/sample_submission.csv')[['customer_ID']]\nsub['customer_ID_hash'] = sub['customer_ID'].str[-16:].str.hex_to_int().astype('int64')\nsub = sub.set_index('customer_ID_hash')\nsub = sub.merge(test[['prediction']], left_index=True, right_index=True, how='left')\nsub = sub.reset_index(drop=True)\n\n# DISPLAY PREDICTIONS\nsub.to_csv(f'submission.csv',index=False)\nprint('Submission file shape is', sub.shape )","metadata":{"_kg_hide-input":true,"id":"xhVpBT6LnSeX","execution":{"iopub.status.busy":"2022-07-20T13:58:12.170410Z","iopub.execute_input":"2022-07-20T13:58:12.170763Z","iopub.status.idle":"2022-07-20T13:58:13.285687Z","shell.execute_reply.started":"2022-07-20T13:58:12.170729Z","shell.execute_reply":"2022-07-20T13:58:13.282208Z"},"trusted":true},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"### 🏁 Conclusion for Hypothesis #4\n\nSubmission with the new features has a score **0.796** on Public LB, and that's more than **0.796** for baseline solution, based on 4th digit after point. Which is not shown ;-).   \nHint - you can check that ranking despite LB rounding to 3 digits in Edit mode -> Competitions. It actually shows what notebook version had maximum Public LB score WITHOUT rounding to 3 digits under the hood (BEST SCORE vs. LATEST SCORE).\n\nWe've got a small improvement both on Local CV and Public LB, as result - for this train-test datasets **we can accept** Hypothesis #4: *Changes in statement dates for 13 month observation period has additional information correlated with default*   ","metadata":{"id":"L5OarQvMnSeX","_kg_hide-input":false}},{"cell_type":"markdown","source":"### 🚀 Useful links with data and feature enrichment guides   \n\n#### [Guide #1 How to improve accuracy of Kaggle TOP1 leaderboard notebook in 10 minutes](https://www.kaggle.com/code/romaupgini/how-to-find-external-data-for-1-private-lb-4-50)\n#### [Guide #2 Zero feature engineering with low-code libraries: Upgini + PyCaret](https://www.kaggle.com/code/romaupgini/zero-feature-engineering-with-upgini-pycaret)\n#### [Guide #3 How to improve accuracy of Multivariate Time Series kernel from external features & data](https://www.kaggle.com/code/romaupgini/guide-external-data-features-for-multivariatets)  \n\n\n#### Happy kaggling! \n<sup>😔 Found error in the library or a bug in notebook code? Our bad! <a href=\"https://github.com/upgini/upgini/issues/new?assignees=&title=readme%2Fbug\">\nPlease report it here.</a></sup>","metadata":{"id":"g9VDS2RMnSeX"}}]}