eda_ml_ea20.ipynb

{
 "cells": [
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "# DATA ANALYSIS MODULE",
   "id": "c8bf8a2b1544f1f4"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-09-23T03:19:12.879689Z",
     "start_time": "2024-09-23T03:19:12.866722Z"
    }
   },
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "import numpy as np\n",
    "from itertools import product\n",
    "\n",
    "# Define the country code at the beginning of the notebook\n",
    "country = 'EA20'  # Change to 'GR', 'DE' OR 'EA20' as needed\n",
    "lower_country = country.lower()\n"
   ],
   "id": "11921dec24b387ef",
   "outputs": [],
   "execution_count": 3
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "# EDA",
   "id": "4438981b61a539be"
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "## ***DATA VISUALIZATION COMPONENT***",
   "id": "bdb07d40fdd517e6"
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Timelines",
   "id": "91f842b49b34acb"
  },
  {
   "cell_type": "code",
   "id": "initial_id",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "# Read the column.md file and store it in a dictionary\n",
    "column_info = {}\n",
    "with open(f'columns_{lower_country}.md', 'r', encoding='utf-8') as file:\n",
    "    for line in file:\n",
    "        line = line.strip()\n",
    "        if line:\n",
    "            col, desc = line.split(': ', 1)  # Use maxsplit=1 to avoid splitting issues\n",
    "            column_info[col] = {'feature': col, 'description': desc}\n",
    "\n",
    "# Read the data from the CSV file\n",
    "data = pd.read_csv(f'INFUSED_DATA/merged_data_{country}.csv')\n",
    "\n",
    "# Convert the TIME_PERIOD column to datetime\n",
    "data['TIME_PERIOD'] = pd.to_datetime(data['TIME_PERIOD'])\n",
    "\n",
    "# Define the number of plots per set and filter only columns that exist in the data\n",
    "plots_per_set = 6\n",
    "filtered_column_info = {k: v for k, v in column_info.items() if k in data.columns}\n",
    "total_features = len(filtered_column_info)\n",
    "\n",
    "# Create the plots\n",
    "for i in range(0, total_features, plots_per_set):\n",
    "    fig, axes = plt.subplots(nrows=3, ncols=2, figsize=(15, 15))\n",
    "    axes = axes.flatten()\n",
    "    \n",
    "    for j, (feature, info) in enumerate(list(filtered_column_info.items())[i:i+plots_per_set]):\n",
    "        ax = axes[j]\n",
    "        ax.plot(data['TIME_PERIOD'], data[feature])\n",
    "        ax.set_title(info['description'])\n",
    "        ax.set_xlabel('TIME_PERIOD')\n",
    "        ax.set_ylabel(info['feature'])\n",
    "    \n",
    "    # Remove empty subplots if there are less than 6 features in this set\n",
    "    for k in range(j+1, plots_per_set):\n",
    "        fig.delaxes(axes[k])\n",
    "    \n",
    "    plt.tight_layout()\n",
    "    plt.show()\n"
   ],
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Distribution diagrams",
   "id": "d91595bc949e1def"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "# Read the column.md file and store it in a dictionary\n",
    "column_info = {}\n",
    "with open(f'columns_{lower_country}.md', 'r', encoding='utf-8') as file:\n",
    "    for line in file:\n",
    "        line = line.strip()\n",
    "        if line:\n",
    "            col, desc = line.split(': ', 1)  # Use maxsplit=1 to avoid splitting issues\n",
    "            column_info[col] = {'feature': col, 'description': desc}\n",
    "\n",
    "# Read the data from the CSV file\n",
    "data = pd.read_csv(f'INFUSED_DATA/merged_data_{country}.csv')\n",
    "\n",
    "# Convert the TIME_PERIOD column to datetime\n",
    "data['TIME_PERIOD'] = pd.to_datetime(data['TIME_PERIOD'])\n",
    "\n",
    "# Define the number of plots per set and filter only columns that exist in the data\n",
    "plots_per_set = 6\n",
    "filtered_column_info = {k: v for k, v in column_info.items() if k in data.columns}\n",
    "total_features = len(filtered_column_info)\n",
    "\n",
    "# Create the plots\n",
    "for i in range(0, total_features, plots_per_set):\n",
    "    fig, axes = plt.subplots(nrows=3, ncols=2, figsize=(15, 15))\n",
    "    axes = axes.flatten()\n",
    "    \n",
    "    for j, (feature, info) in enumerate(list(filtered_column_info.items())[i:i+plots_per_set]):\n",
    "        ax = axes[j]\n",
    "        sns.histplot(data[feature], kde=True, ax=ax)  # Use seaborn to create histogram with density line\n",
    "        ax.set_title(info['description'])\n",
    "        ax.set_xlabel(feature)\n",
    "        ax.set_ylabel('Frequency')\n",
    "    \n",
    "    # Remove empty subplots if there are less than 6 features in this set\n",
    "    for k in range(j+1, plots_per_set):\n",
    "        fig.delaxes(axes[k])\n",
    "    \n",
    "    plt.tight_layout()\n",
    "    plt.show()\n"
   ],
   "id": "9b3007db0880e90",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Boxplots",
   "id": "ec0528b8ad28c877"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "# Read the column.md file and store it in a dictionary\n",
    "column_info = {}\n",
    "with open(f'columns_{lower_country}.md', 'r', encoding='utf-8') as file:\n",
    "    for line in file:\n",
    "        line = line.strip()\n",
    "        if line:\n",
    "            col, desc = line.split(': ', 1)  # Use maxsplit=1 to avoid splitting issues\n",
    "            column_info[col] = {'feature': col, 'description': desc}\n",
    "\n",
    "# Read the data from the CSV file\n",
    "data = pd.read_csv(f'INFUSED_DATA/merged_data_{country}.csv')\n",
    "\n",
    "# Convert the TIME_PERIOD column to datetime\n",
    "data['TIME_PERIOD'] = pd.to_datetime(data['TIME_PERIOD'])\n",
    "\n",
    "# Function to calculate outliers percentage\n",
    "def calculate_outliers(data, feature, threshold=1.5):\n",
    "    q1 = data[feature].quantile(0.25)\n",
    "    q3 = data[feature].quantile(0.75)\n",
    "    iqr = q3 - q1\n",
    "    lower_bound = q1 - threshold * iqr\n",
    "    upper_bound = q3 + threshold * iqr\n",
    "\n",
    "    outliers = data[(data[feature] < lower_bound) | (data[feature] > upper_bound)]\n",
    "    outliers_count = outliers.shape[0]\n",
    "    total_count = data.shape[0]\n",
    "    outliers_percentage = (outliers_count / total_count) * 100\n",
    "    \n",
    "    return outliers_count, outliers_percentage\n",
    "\n",
    "# Function to plot outliers\n",
    "def plot_outliers(data, features, title):\n",
    "    plots_per_page = 6\n",
    "    rows, cols = 3, 2\n",
    "\n",
    "    for i in range(0, len(features), plots_per_page):\n",
    "        fig, axs = plt.subplots(rows, cols, figsize=(15, 15))\n",
    "        axs = axs.flatten()\n",
    "\n",
    "        for j, feature in enumerate(features[i:i+plots_per_page]):\n",
    "            if feature in data.columns:  # Check if the feature exists in the data\n",
    "                ax = axs[j]\n",
    "                sns.boxplot(x=data[feature], ax=ax)\n",
    "                outliers_count, outliers_percentage = calculate_outliers(data, feature)\n",
    "                ax.set_title(f'Boxplot for {feature}\\nOutliers: {outliers_count} ({outliers_percentage:.2f}%)')\n",
    "                ax.set_xlabel(feature)\n",
    "\n",
    "        # Remove any unused subplots\n",
    "        for k in range(j+1, len(axs)):\n",
    "            fig.delaxes(axs[k])\n",
    "\n",
    "        plt.suptitle(title, fontsize=16)\n",
    "        plt.tight_layout(rect=[0, 0, 1, 0.96])  # Adjust layout to make space for the main title\n",
    "        plt.show()\n",
    "        plt.close(fig)  # Close the figure to free memory\n",
    "\n",
    "# List of features for outliers analysis\n",
    "all_features = list(column_info.keys())\n",
    "\n",
    "# Filter features to include only those present in the data\n",
    "filtered_features = [feature for feature in all_features if feature in data.columns]\n",
    "\n",
    "# Plot outliers\n",
    "plot_outliers(data, filtered_features, 'Outliers Detection for All Features')\n"
   ],
   "id": "8ad67c7e117d47bb",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Data Transformation",
   "id": "222d1bb9dede6620"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-24T05:20:34.475122Z",
     "start_time": "2024-08-24T05:20:34.337202Z"
    }
   },
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "# Read the data from the updated CSV file\n",
    "data = pd.read_csv(f'INFUSED_DATA/merged_data_{country}.csv')\n",
    "\n",
    "# List of trend and VideoCount variables to convert to binary\n",
    "variables = [\n",
    "    'Russia_Ukraine_War_Trend', \n",
    "    'Israel_Gaza_War_Trend', \n",
    "    'Financial_Crisis_Trend', \n",
    "    'Climate_Change_Trend', \n",
    "    'Covid19_Trend',\n",
    "    'Russia_Ukraine_War_VideoCount', \n",
    "    'Israel_Gaza_War_VideoCount', \n",
    "    'Financial_Crisis_VideoCount', \n",
    "    'Climate_Change_VideoCount', \n",
    "    'Covid19_VideoCount'\n",
    "]\n",
    "\n",
    "# Function to convert to binary based on mean\n",
    "def convert_to_binary(df, col):\n",
    "    threshold = df[col].mean()\n",
    "    binary_col = f'{col}_Binary'\n",
    "    df[binary_col] = np.where(df[col] >= threshold, 1, 0)\n",
    "\n",
    "# Apply the function to each variable\n",
    "for var in variables:\n",
    "    convert_to_binary(data, var)\n",
    "\n",
    "# Save the updated data to a new CSV file\n",
    "data.to_csv(f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv', index=False)\n"
   ],
   "id": "e2e48addb67cea7f",
   "outputs": [],
   "execution_count": 5
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Correlation Analysis",
   "id": "bff2299461d3bcc4"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns\n",
    "\n",
    "# Read the data from the updated CSV file\n",
    "data = pd.read_csv(f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv')\n",
    "\n",
    "# Convert the TIME_PERIOD column to datetime\n",
    "data['TIME_PERIOD'] = pd.to_datetime(data['TIME_PERIOD'])\n",
    "\n",
    "# Compute and plot correlation matrix for original data\n",
    "original_data_no_time_period = data.drop(columns=['TIME_PERIOD'])\n",
    "correlation_matrix_original = original_data_no_time_period.corr()\n",
    "plt.figure(figsize=(20, 16))\n",
    "sns.heatmap(correlation_matrix_original, annot=True, cmap='coolwarm', linewidths=.5)\n",
    "plt.title('Heatmap of Feature Correlations')\n",
    "plt.show()\n",
    "\n",
    "# Extracting strong correlations and saving to CSV\n",
    "def extract_strong_correlations(corr_matrix, threshold=0.60):\n",
    "    corr_pairs = []\n",
    "    for i in range(len(corr_matrix.columns)):\n",
    "        for j in range(i+1, len(corr_matrix.columns)):\n",
    "            if abs(corr_matrix.iloc[i, j]) >= threshold:\n",
    "                corr_pairs.append((corr_matrix.columns[i], corr_matrix.columns[j], corr_matrix.iloc[i, j]))\n",
    "    return corr_pairs\n",
    "\n",
    "# Extract strong correlations\n",
    "strong_correlations = extract_strong_correlations(correlation_matrix_original, threshold=0.60)\n",
    "\n",
    "# Create a DataFrame for strong correlations\n",
    "strong_correlations_df = pd.DataFrame(strong_correlations, columns=['Feature 1', 'Feature 2', 'Correlation'])\n",
    "\n",
    "# Split the DataFrame into positive and negative correlations\n",
    "positive_corr_df = strong_correlations_df[strong_correlations_df['Correlation'] > 0].copy()\n",
    "negative_corr_df = strong_correlations_df[strong_correlations_df['Correlation'] < 0].copy()\n",
    "\n",
    "# Sort by absolute value of correlation\n",
    "positive_corr_df = positive_corr_df.reindex(positive_corr_df['Correlation'].abs().sort_values(ascending=False).index)\n",
    "negative_corr_df = negative_corr_df.reindex(negative_corr_df['Correlation'].abs().sort_values(ascending=False).index)\n",
    "\n",
    "# Concatenate the positive and negative correlations\n",
    "sorted_strong_correlations_df = pd.concat([positive_corr_df, negative_corr_df])\n",
    "\n",
    "# Save the strong correlations to a CSV file\n",
    "sorted_strong_correlations_df.to_csv(f'EDA_ML_DATA/{country}/strong_correlations.csv', index=False)\n"
   ],
   "id": "4c109bae4d01eb53",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Scatter Plots",
   "id": "8db764f1fb6c7cf9"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "import pandas as pd\n",
    "\n",
    "# Load the strong correlations DataFrame\n",
    "strong_corr_features = pd.read_csv(f'EDA_ML_DATA/{country}/strong_correlations.csv')\n",
    "\n",
    "# Define the exclusion keywords\n",
    "exclude_keywords = ['BS', 'CISS', 'EU', 'ESTR', 'UN', 'ECB', 'Close', 'Volume', 'Deaths', 'Cases', 'Trend', 'VideoCount', 'Binary']\n",
    "\n",
    "# Function to check if both features contain the same keyword\n",
    "def contains_same_keyword(feature1, feature2, keywords):\n",
    "    for keyword in keywords:\n",
    "        if keyword in feature1 and keyword in feature2:\n",
    "            return True\n",
    "    return False\n",
    "\n",
    "# Filter out rows where both 'Feature 1' and 'Feature 2' contain the same keyword\n",
    "strong_corr_features = strong_corr_features[\n",
    "    ~strong_corr_features.apply(lambda row: contains_same_keyword(row['Feature 1'], row['Feature 2'], exclude_keywords), axis=1)\n",
    "]\n",
    "\n",
    "# Define strong correlations (correlation > 0.75 or correlation < -0.75)\n",
    "strong_corr_features = strong_corr_features[(strong_corr_features['Correlation'] > 0.75) | (strong_corr_features['Correlation'] < -0.75)]\n",
    "\n",
    "# Split the features into positive and negative correlations\n",
    "positive_corr_features = strong_corr_features[strong_corr_features['Correlation'] >= 0.75].sort_values(by='Correlation', ascending=False)\n",
    "negative_corr_features = strong_corr_features[strong_corr_features['Correlation'] <= -0.75].sort_values(by='Correlation')\n",
    "\n",
    "# Function to plot scatter plots with regression lines\n",
    "def plot_correlations(corr_features, title_prefix, data):\n",
    "    plots_per_page = 6\n",
    "    rows, cols = 3, 2\n",
    "\n",
    "    for i in range(0, len(corr_features), plots_per_page):\n",
    "        fig, axs = plt.subplots(rows, cols, figsize=(15, 15))\n",
    "        axs = axs.flatten()\n",
    "\n",
    "        for j, (index, row) in enumerate(corr_features.iloc[i:i+plots_per_page].iterrows()):\n",
    "            ax = axs[j]\n",
    "            sns.scatterplot(data=data, x=row['Feature 1'], y=row['Feature 2'], ax=ax)\n",
    "            sns.regplot(data=data, x=row['Feature 1'], y=row['Feature 2'], scatter=False, color='red', ax=ax)\n",
    "            ax.set_title(f'{title_prefix}: {row[\"Feature 1\"]} vs {row[\"Feature 2\"]} (Correlation: {row[\"Correlation\"]:.2f})')\n",
    "            ax.set_xlabel(row['Feature 1'])\n",
    "            ax.set_ylabel(row['Feature 2'])\n",
    "\n",
    "        # Remove any unused subplots\n",
    "        for j in range(len(corr_features.iloc[i:i+plots_per_page]), len(axs)):\n",
    "            fig.delaxes(axs[j])\n",
    "\n",
    "        plt.tight_layout()\n",
    "        plt.show()\n",
    "        plt.close(fig)  # Close the figure to free memory\n",
    "\n",
    "# Load the data\n",
    "data = pd.read_csv(f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv')\n",
    "\n",
    "# Plot positive correlations\n",
    "plot_correlations(positive_corr_features, 'Positive Correlation', data)\n",
    "\n",
    "# Plot negative correlations\n",
    "plot_correlations(negative_corr_features, 'Negative Correlation', data)\n"
   ],
   "id": "19fd4c79d79919d1",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### Groupwise Correlation Analysis",
   "id": "b1c5de90b2084e8e"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "import seaborn as sns\n",
    "import matplotlib.pyplot as plt\n",
    "from itertools import combinations\n",
    "\n",
    "# Load the data\n",
    "data = pd.read_csv(f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv')\n",
    "\n",
    "# Define feature groups\n",
    "economic_features = [\n",
    "    f'GDP_GROWTH_{country}', f'INFLATION_{country}_VALUE', f'GOVDEBT_{country}',\n",
    "    'ECB_EXCHANGE_RATES_USD_EUR', 'ECB_EXCHANGE_RATES_CNY_EUR',\n",
    "    'ESTR_EU000A2QQF08_CI', 'ESTR_EU000A2X2A25_NT', 'ESTR_EU000A2X2A25_TT',\n",
    "    f'CISS_{country}_SS_CIN', 'CISS_EA20_SS_BM', 'CISS_EA20_SS_FI',\n",
    "    'CISS_EA20_SS_FX', 'CISS_EA20_SS_MM',\n",
    "    f'BS_{country}_CSMCI', f'BS_{country}_FS_NY', f'BS_{country}_GES_NY', f'BS_{country}_MP_NY',\n",
    "    f'BS_{country}_PT_NY', f'BS_{country}_SV_NY', f'UN_{country}_T_TOT', f'UN_{country}_M_TOT', \n",
    "    f'UN_{country}_F_TOT', f'UN_{country}_GT25_TOT', f'UN_{country}_LE25_TOT', 'Oil_Close_Price',\n",
    "    'Oil_Volume', 'Gas_Close_Price', 'Gas_Volume', 'Corn_Close_Price', 'Corn_Volume',\n",
    "    'Wheat_Close_Price', 'Wheat_Volume'\n",
    "]\n",
    "covid_features = [col for col in data.columns if 'Covid_' in col]\n",
    "trend_features = [col for col in data.columns if col.endswith('Trend') or col.endswith('Trend_Binary')]\n",
    "videocount_features = [col for col in data.columns if col.endswith('VideoCount') or col.endswith('VideoCount_Binary')]\n",
    "\n",
    "# Function to create heatmap for a set of features\n",
    "def plot_heatmap(features, title):\n",
    "    subset = data[features].dropna()\n",
    "    correlation_matrix = subset.corr()\n",
    "    plt.figure(figsize=(20, 10))  # Adjust the figure size as needed\n",
    "    sns.heatmap(correlation_matrix, annot=True, cmap='coolwarm', linewidths=.5)\n",
    "    plt.title(f'Heatmap of {title}')\n",
    "    plt.show()\n",
    "\n",
    "# Plot heatmap for economic features separately\n",
    "plot_heatmap(economic_features, 'Economic Features')\n",
    "\n",
    "# Plot heatmap for covid features separately\n",
    "plot_heatmap(covid_features, 'Covid Features')\n",
    "\n",
    "# Plot heatmap for trend features separately\n",
    "plot_heatmap(trend_features, 'Trend Features')\n",
    "\n",
    "# Plot heatmap for videocount features separately\n",
    "plot_heatmap(videocount_features, 'VideoCount Features')\n",
    "\n",
    "# Function to create heatmap of correlations between two feature groups\n",
    "def plot_intergroup_heatmap(group1, group2, title):\n",
    "    subset = data[group1 + group2].dropna()\n",
    "    correlation_matrix = subset.corr().loc[group1, group2]\n",
    "    plt.figure(figsize=(20, 10))  # Adjust the figure size as needed\n",
    "    sns.heatmap(correlation_matrix, annot=True, cmap='coolwarm', linewidths=.5)\n",
    "    plt.title(f'Heatmap of {title}')\n",
    "    plt.show()\n",
    "\n",
    "# Define feature groups\n",
    "groups = {\n",
    "    'Economic Features': economic_features,\n",
    "    'Covid Features': covid_features,\n",
    "    'Trend Features': trend_features,\n",
    "    'VideoCount Features': videocount_features\n",
    "}\n",
    "\n",
    "# Create heatmaps for correlations between different feature groups\n",
    "group_pairs = list(combinations(groups.items(), 2))\n",
    "\n",
    "for (group1_name, group1_features), (group2_name, group2_features) in group_pairs:\n",
    "    plot_intergroup_heatmap(group1_features, group2_features, f'{group1_name} vs {group2_name}')\n"
   ],
   "id": "43c4fa456ff8eb65",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "## ***PREPARATION AND SELECTION OF VARIABLE GROUPS FOR THE MACHINE LEARNING STAGE***",
   "id": "ebaac1a69cd6d776"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "from statsmodels.stats.outliers_influence import variance_inflation_factor\n",
    "from statsmodels.tools.tools import add_constant\n",
    "import ast\n",
    "\n",
    "# Function to calculate VIF\n",
    "def calculate_vif(df, features):\n",
    "    X = df[features]\n",
    "\n",
    "    # Ensure columns are numeric\n",
    "    for col in features:\n",
    "        if not pd.api.types.is_numeric_dtype(X[col]):\n",
    "            raise ValueError(f\"The column {col} is not numeric.\")\n",
    "\n",
    "    X = X.dropna()  # Remove rows with NaN values\n",
    "    X = add_constant(X)\n",
    "\n",
    "    vif_data = pd.DataFrame()\n",
    "    vif_data[\"Feature\"] = X.columns\n",
    "    vif_data[\"VIF\"] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]\n",
    "    return vif_data\n",
    "\n",
    "# Read the CSV file with correlations\n",
    "file_path = f'EDA_ML_DATA/{country}/strong_correlations.csv'\n",
    "df = pd.read_csv(file_path)\n",
    "\n",
    "# List of all unique variables from 'Feature 1' column\n",
    "unique_variables = df['Feature 1'].unique()\n",
    "\n",
    "# Create all possible groups of dependent-independent variables\n",
    "groups = []\n",
    "\n",
    "for dependent_var in unique_variables:\n",
    "    # Select rows where 'Feature 1' is the dependent variable and the correlation is above 0.60\n",
    "    relevant_rows = df[(df['Feature 1'] == dependent_var) & (df['Correlation'].abs() > 0.60)]\n",
    "\n",
    "    # Get independent variables\n",
    "    independent_vars = relevant_rows['Feature 2'].tolist()\n",
    "\n",
    "    if independent_vars:\n",
    "        groups.append((dependent_var, independent_vars))\n",
    "\n",
    "# Convert groups to DataFrame\n",
    "groups_df = pd.DataFrame(groups, columns=['Dependent Variable', 'Independent Variables'])\n",
    "\n",
    "# Save the result to a CSV file\n",
    "output_file_path = f'EDA_ML_DATA/{country}/dependent_independent_groups.csv'\n",
    "groups_df.to_csv(output_file_path, index=False)\n",
    "\n",
    "print(f'The dependent-independent groups have been saved to {output_file_path}')\n",
    "\n",
    "# Read the CSV file with the groups\n",
    "groups_file_path = f'EDA_ML_DATA/{country}/dependent_independent_groups.csv'\n",
    "groups_df = pd.read_csv(groups_file_path)\n",
    "\n",
    "# Read the CSV file with correlations\n",
    "correlations_file_path = f'EDA_ML_DATA/{country}/strong_correlations.csv'\n",
    "correlations_df = pd.read_csv(correlations_file_path)\n",
    "\n",
    "# Read the dataset with actual values\n",
    "dataset_file_path = f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv'\n",
    "full_df = pd.read_csv(dataset_file_path)\n",
    "\n",
    "# Process to check the independent variables of each group\n",
    "results = []\n",
    "\n",
    "for _, row in groups_df.iterrows():\n",
    "    dependent_var = row['Dependent Variable']\n",
    "    independent_vars = ast.literal_eval(row['Independent Variables'])  # Convert string list to list\n",
    "\n",
    "    # Check correlations\n",
    "    correlated_pairs = []\n",
    "    for i, var1 in enumerate(independent_vars):\n",
    "        for var2 in independent_vars[i+1:]:\n",
    "            corr_row = correlations_df[\n",
    "                ((correlations_df['Feature 1'] == var1) & (correlations_df['Feature 2'] == var2)) |\n",
    "                ((correlations_df['Feature 1'] == var2) & (correlations_df['Feature 2'] == var1))\n",
    "            ]\n",
    "            if not corr_row.empty and abs(corr_row['Correlation'].values[0]) > 0.75:\n",
    "                correlated_pairs.append((var1, var2, corr_row['Correlation'].values[0]))\n",
    "\n",
    "    # Calculate VIF using the actual values from full_df\n",
    "    try:\n",
    "        vif_result = calculate_vif(full_df, independent_vars)\n",
    "    except Exception as e:\n",
    "        vif_result = pd.DataFrame({'Feature': independent_vars, 'VIF': [str(e)]*len(independent_vars)})\n",
    "\n",
    "    # Save results\n",
    "    results.append({\n",
    "        'Dependent Variable': dependent_var,\n",
    "        'Independent Variables': independent_vars,\n",
    "        'Correlated Pairs': correlated_pairs,\n",
    "        'VIF Result': vif_result\n",
    "    })\n",
    "\n",
    "# Create final DataFrame for all results\n",
    "final_results = []\n",
    "\n",
    "for result in results:\n",
    "    dependent_var = result['Dependent Variable']\n",
    "    for i, row in result['VIF Result'].iterrows():\n",
    "        feature = row['Feature']\n",
    "        vif = row['VIF']\n",
    "        correlated_with = [pair for pair in result['Correlated Pairs'] if feature in pair[:2]]\n",
    "        final_results.append({\n",
    "            'Dependent Variable': dependent_var,\n",
    "            'Feature': feature,\n",
    "            'VIF': vif,\n",
    "            'Correlated With': correlated_with\n",
    "        })\n",
    "\n",
    "final_results_df = pd.DataFrame(final_results)\n",
    "\n",
    "# Save the final file\n",
    "final_output_file = f'EDA_ML_DATA/{country}/final_vif_correlations_results.csv'\n",
    "final_results_df.to_csv(final_output_file, index=False)\n",
    "\n",
    "print(f\"The final results have been saved to {final_output_file}\")\n"
   ],
   "id": "42694a7190b3a458",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "##### VIF-High Correlation Testing",
   "id": "3347982b5f09d96b"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "from statsmodels.stats.outliers_influence import variance_inflation_factor\n",
    "from statsmodels.tools.tools import add_constant\n",
    "\n",
    "# Function to calculate VIF\n",
    "def calculate_vif(df, features):\n",
    "    X = df[features]\n",
    "\n",
    "    # Ensure columns are numeric\n",
    "    for col in features:\n",
    "        if not pd.api.types.is_numeric_dtype(X[col]):\n",
    "            raise ValueError(f\"The column {col} is not numeric.\")\n",
    "\n",
    "    X = X.dropna()  # Remove rows with NaN values\n",
    "    X = add_constant(X)\n",
    "\n",
    "    vif_data = pd.DataFrame()\n",
    "    vif_data[\"Feature\"] = X.columns\n",
    "    vif_data[\"VIF\"] = [variance_inflation_factor(X.values, i) for i in range(X.shape[1])]\n",
    "    return vif_data\n",
    "\n",
    "# Define dependent and independent variables for Team-1\n",
    "#dependent_var = f'UN_{country}_T_TOT'  # Replace with your dependent variable\n",
    "#independent_vars = ['Covid19_Trend_Binary', 'ECB_EXCHANGE_RATES_USD_EUR', f'Covid_{country}_Cumulative_Deaths', 'Oil_Close_Price', 'Corn_Close_Price']  # Replace with your independent variables\n",
    "\n",
    "# Define dependent and independent variables for Team-2\n",
    "dependent_var = f'GOVDEBT_{country}'  # Replace with your dependent variable Israel_Gaza_War_VideoCount_Binary\n",
    "independent_vars = ['BS_EA20_MP_NY', 'UN_EA20_T_TOT', f'Covid_{country}_Cumulative_Deaths', 'Oil_Close_Price']   # Replace with your independent variables\n",
    "\n",
    "# Define dependent and independent variables for Team-3\n",
    "#dependent_var = f'INFLATION_{country}_VALUE'  # Replace with your dependent variable\n",
    "#independent_vars = ['Russia_Ukraine_War_VideoCount_Binary', 'Oil_Close_Price', f'Covid_{country}_Cumulative_Deaths', 'ECB_EXCHANGE_RATES_CNY_EUR', f'BS_{country}_FS_NY', f'UN_{country}_T_TOT', 'Covid19_Trend_Binary', 'ECB_EXCHANGE_RATES_USD_EUR']  # Replace with your independent variables\n",
    "\n",
    "# Read the dataset with actual values\n",
    "dataset_file_path = f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv'\n",
    "full_df = pd.read_csv(dataset_file_path)\n",
    "\n",
    "# Check correlations between independent variables\n",
    "correlations_df = pd.read_csv(f'EDA_ML_DATA/{country}/strong_correlations.csv')\n",
    "\n",
    "correlated_pairs = []\n",
    "for i, var1 in enumerate(independent_vars):\n",
    "    for var2 in independent_vars[i+1:]:\n",
    "        corr_row = correlations_df[\n",
    "            ((correlations_df['Feature 1'] == var1) & (correlations_df['Feature 2'] == var2)) |\n",
    "            ((correlations_df['Feature 1'] == var2) & (correlations_df['Feature 2'] == var1))\n",
    "        ]\n",
    "        if not corr_row.empty and abs(corr_row['Correlation'].values[0]) > 0.75:\n",
    "            correlated_pairs.append((var1, var2, corr_row['Correlation'].values[0]))\n",
    "\n",
    "# Calculate VIF using the actual values from full_df\n",
    "try:\n",
    "    vif_result = calculate_vif(full_df, independent_vars)\n",
    "except Exception as e:\n",
    "    vif_result = pd.DataFrame({'Feature': independent_vars, 'VIF': [str(e)]*len(independent_vars)})\n",
    "\n",
    "# Save results\n",
    "results = [{\n",
    "    'Dependent Variable': dependent_var,\n",
    "    'Independent Variables': independent_vars,\n",
    "    'Correlated Pairs': correlated_pairs,\n",
    "    'VIF Result': vif_result\n",
    "}]\n",
    "\n",
    "# Create final DataFrame for the results\n",
    "final_results = []\n",
    "\n",
    "for result in results:\n",
    "    dependent_var = result['Dependent Variable']\n",
    "    for i, row in result['VIF Result'].iterrows():\n",
    "        feature = row['Feature']\n",
    "        vif = row['VIF']\n",
    "        correlated_with = [pair for pair in result['Correlated Pairs'] if feature in pair[:2]]\n",
    "        final_results.append({\n",
    "            'Dependent Variable': dependent_var,\n",
    "            'Feature': feature,\n",
    "            'VIF': vif,\n",
    "            'Correlated With': correlated_with\n",
    "        })\n",
    "\n",
    "final_results_df = pd.DataFrame(final_results)\n",
    "\n",
    "# Save the final file\n",
    "final_output_file = f'EDA_ML_DATA/{country}/testing_vif_correlations_results.csv'\n",
    "final_results_df.to_csv(final_output_file, index=False)\n",
    "\n",
    "print(f\"The testing results have been saved to {final_output_file}\")\n"
   ],
   "id": "27a981c2f0372285",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "# HYPOTHESIS TESTING",
   "id": "4180233752f5c89b"
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "## ***Linear Regression - Random Forests - Gradient Boosting***",
   "id": "2aa7e82ee586af3e"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "from sklearn.linear_model import LinearRegression\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.metrics import mean_absolute_error, mean_squared_error\n",
    "import numpy as np\n",
    "import statsmodels.api as sm\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.ensemble import RandomForestRegressor, GradientBoostingRegressor\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.model_selection import cross_val_score\n",
    "from sklearn.metrics import make_scorer\n",
    "\n",
    "# Read the dataset with actual values\n",
    "dataset_file_path = f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv'\n",
    "df = pd.read_csv(dataset_file_path)\n",
    "\n",
    "# Convert TIME_PERIOD column to datetime\n",
    "df['TIME_PERIOD'] = pd.to_datetime(df['TIME_PERIOD'])\n",
    "\n",
    "# Define dependent and independent variables for Team-1  \n",
    "dependent_var = f'UN_{country}_T_TOT'  # Replace with your dependent variable \n",
    "independent_vars = ['Covid19_Trend_Binary', 'ECB_EXCHANGE_RATES_USD_EUR', 'Oil_Close_Price', f'Covid_{country}_Cumulative_Deaths',  'Corn_Close_Price']  # Replace with your independent variables\n",
    "\n",
    "# Define dependent and independent variables for Team-2\n",
    "#  \n",
    "#dependent_var = f'GOVDEBT_{country}'  # Replace with your dependent variable\n",
    "#independent_vars = ['BS_EA20_MP_NY', f'Covid_{country}_Cumulative_Deaths', 'UN_EA20_T_TOT', 'Oil_Close_Price']   # Replace with your independent variables\n",
    "\n",
    "# Create X and y\n",
    "X = df[independent_vars]\n",
    "y = df[dependent_var]\n",
    "\n",
    "# Split the data into training and testing sets using random order\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)\n",
    "\n",
    "# Z-score normalization for all models\n",
    "scaler_X = StandardScaler()\n",
    "X_train_scaled = scaler_X.fit_transform(X_train)\n",
    "X_test_scaled = scaler_X.transform(X_test)\n",
    "\n",
    "scaler_y = StandardScaler()\n",
    "y_train_scaled = scaler_y.fit_transform(y_train.values.reshape(-1, 1))\n",
    "y_test_scaled = scaler_y.transform(y_test.values.reshape(-1, 1))\n",
    "\n",
    "# Linear Regression Model\n",
    "model_lr = LinearRegression()\n",
    "model_lr.fit(X_train_scaled, y_train_scaled)\n",
    "\n",
    "# Fit the model with statsmodels for detailed analysis\n",
    "X_train_scaled_sm = sm.add_constant(X_train_scaled)  # Add constant term\n",
    "model_sm = sm.OLS(y_train_scaled, X_train_scaled_sm).fit()\n",
    "\n",
    "# Print coefficients, p-values, and check for significance\n",
    "print(\"Linear Regression Coefficients and p-values:\")\n",
    "coefficients = model_sm.params\n",
    "p_values = model_sm.pvalues\n",
    "for coef, p_val, var in zip(coefficients, p_values, ['const'] + independent_vars):\n",
    "    significance = \"Significant\" if p_val < 0.05 else \"Not Significant\"\n",
    "    print(f\"{var}: Coefficient = {coef:.4f}, p-value = {p_val:.4f} ({significance})\")\n",
    "\n",
    "# Make predictions with Linear Regression\n",
    "predictions_scaled_lr = model_lr.predict(X_test_scaled)\n",
    "predictions_lr = scaler_y.inverse_transform(predictions_scaled_lr)\n",
    "\n",
    "# Calculate evaluation metrics for Linear Regression\n",
    "r_squared_lr = model_lr.score(X_test_scaled, y_test_scaled)\n",
    "mae_lr = mean_absolute_error(y_test, predictions_lr)\n",
    "mse_lr = mean_squared_error(y_test, predictions_lr)\n",
    "rmse_lr = np.sqrt(mse_lr)\n",
    "\n",
    "# Calculate mean of actual values\n",
    "mean_y_test = np.mean(y_test)\n",
    "\n",
    "# Calculate percentage metrics\n",
    "mae_lr_percent = (mae_lr / mean_y_test) * 100\n",
    "rmse_lr_percent = (rmse_lr / mean_y_test) * 100\n",
    "mse_lr_percent = (mse_lr / (mean_y_test ** 2)) * 100\n",
    "\n",
    "print(f\"Linear Regression - R-squared: {r_squared_lr}\")\n",
    "print(f\"Linear Regression - Mean Absolute Error (MAE): {mae_lr} ({mae_lr_percent:.2f}%)\")\n",
    "print(f\"Linear Regression - Mean Squared Error (MSE): {mse_lr} ({mse_lr_percent:.2f}%)\")\n",
    "print(f\"Linear Regression - Root Mean Squared Error (RMSE): {rmse_lr} ({rmse_lr_percent:.2f}%)\")\n",
    "mean_y_test = np.mean(y_test)\n",
    "print(f\"Mean of dependent variable (y_test): {mean_y_test:.4f}\")\n",
    "\n",
    "# Create a dataframe to hold the actual and predicted values\n",
    "results_df = df[['TIME_PERIOD', dependent_var]].copy()\n",
    "results_df['Predicted_LR'] = np.nan\n",
    "results_df.loc[X_test.index, 'Predicted_LR'] = predictions_lr\n",
    "\n",
    "# Plot the actual values and predicted values over time for Linear Regression\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.plot(results_df['TIME_PERIOD'], results_df[dependent_var], label='Actual Values', color='b', linewidth = 0.8)\n",
    "plt.plot(results_df['TIME_PERIOD'], results_df['Predicted_LR'], label='Predicted Values', color='r', linestyle = '--', marker = 'o', markersize = 1)\n",
    "plt.xlabel('Time')\n",
    "plt.ylabel(dependent_var)\n",
    "plt.title(f'Actual vs Predicted Values for {dependent_var} - Linear Regression')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "# Plotting the importance of variables based on p-values\n",
    "plt.figure(figsize=(10, 6))\n",
    "p_values_dict = {var: p_val for var, p_val in zip(['const'] + independent_vars, p_values)}\n",
    "sorted_p_values = dict(sorted(p_values_dict.items(), key=lambda item: item[1]))\n",
    "\n",
    "bars = plt.bar(sorted_p_values.keys(), -np.log10(list(sorted_p_values.values())), color=['gray' if var == 'const' else 'blue' for var in sorted_p_values.keys()])\n",
    "plt.axhline(y=-np.log10(0.05), color='red', linestyle='--', label='Significance Level (0.05)')\n",
    "plt.xlabel('Variables')\n",
    "plt.ylabel('-log10(p-value)')\n",
    "plt.title('Variable Significance Based on Hypothesis Testing')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.xticks(rotation=45)\n",
    "plt.show()\n",
    "\n",
    "# Random Forest Model using scaled data\n",
    "print(\"\\nRandom Forest Model Results:\")\n",
    "rf_model = RandomForestRegressor(n_estimators=100, random_state=42)\n",
    "rf_model.fit(X_train_scaled, y_train_scaled.ravel())\n",
    "rf_predictions_scaled = rf_model.predict(X_test_scaled)\n",
    "rf_predictions = scaler_y.inverse_transform(rf_predictions_scaled.reshape(-1, 1))\n",
    "\n",
    "# Calculate evaluation metrics for Random Forest\n",
    "rf_mae = mean_absolute_error(y_test, rf_predictions)\n",
    "rf_mse = mean_squared_error(y_test, rf_predictions)\n",
    "rf_rmse = np.sqrt(rf_mse)\n",
    "\n",
    "# Calculate percentage metrics for Random Forest\n",
    "rf_mae_percent = (rf_mae / mean_y_test) * 100\n",
    "rf_rmse_percent = (rf_rmse / mean_y_test) * 100\n",
    "rf_mse_percent = (rf_mse / (mean_y_test ** 2)) * 100\n",
    "\n",
    "print(f\"Random Forest - Mean Absolute Error (MAE): {rf_mae} ({rf_mae_percent:.2f}%)\")\n",
    "print(f\"Random Forest - Mean Squared Error (MSE): {rf_mse} ({rf_mse_percent:.2f}%)\")\n",
    "print(f\"Random Forest - Root Mean Squared Error (RMSE): {rf_rmse} ({rf_rmse_percent:.2f}%)\")\n",
    "\n",
    "# Add RF predictions to the results dataframe\n",
    "results_df['Predicted_RF'] = np.nan\n",
    "results_df.loc[X_test.index, 'Predicted_RF'] = rf_predictions\n",
    "\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.plot(results_df['TIME_PERIOD'], y, label='Actual Values', color='b', linewidth = 0.8)\n",
    "plt.plot(results_df['TIME_PERIOD'], results_df['Predicted_RF'], label='RF Predicted Values', color = 'r', linestyle = '--', marker = 'o', markersize = 1)\n",
    "plt.xlabel('Time')\n",
    "plt.ylabel(dependent_var)\n",
    "plt.title(f'Actual vs Random Forest Predicted Values for {dependent_var}')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "importances_rf = rf_model.feature_importances_\n",
    "indices_rf = np.argsort(importances_rf)[::-1]\n",
    "features_rf = [independent_vars[i] for i in indices_rf]\n",
    "\n",
    "# Print feature importances for Random Forest\n",
    "print(\"\\nFeature Importances based on Gini Index (Random Forest):\")\n",
    "for feature_name, importance in zip(features_rf, importances_rf[indices_rf]):\n",
    "    print(f\"{feature_name}: {importance:.4f}\")\n",
    "\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.bar(features_rf, importances_rf[indices_rf], color='blue')\n",
    "plt.axhline(y=0.05, color='red', linestyle='--', label='Significance Level (0.05)')\n",
    "plt.xlabel('Variables')\n",
    "plt.ylabel('Importance')\n",
    "plt.title('Variable Importance Based on Random Forest')\n",
    "plt.xticks(rotation=45)\n",
    "plt.legend()  # Add legend to show the significance level\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "# Gradient Boosting Model using scaled data\n",
    "print(\"\\nGradient Boosting Model Results:\")\n",
    "gb_model = GradientBoostingRegressor(n_estimators=100, random_state=42)\n",
    "gb_model.fit(X_train_scaled, y_train_scaled.ravel())\n",
    "gb_predictions_scaled = gb_model.predict(X_test_scaled)\n",
    "gb_predictions = scaler_y.inverse_transform(gb_predictions_scaled.reshape(-1, 1))\n",
    "\n",
    "# Calculate evaluation metrics for Gradient Boosting\n",
    "gb_mae = mean_absolute_error(y_test, gb_predictions)\n",
    "gb_mse = mean_squared_error(y_test, gb_predictions)\n",
    "gb_rmse = np.sqrt(gb_mse)\n",
    "\n",
    "# Calculate percentage metrics for Gradient Boosting\n",
    "gb_mae_percent = (gb_mae / mean_y_test) * 100\n",
    "gb_rmse_percent = (gb_rmse / mean_y_test) * 100\n",
    "gb_mse_percent = (gb_mse / (mean_y_test ** 2)) * 100\n",
    "\n",
    "print(f\"Gradient Boosting - Mean Absolute Error (MAE): {gb_mae} ({gb_mae_percent:.2f}%)\")\n",
    "print(f\"Gradient Boosting - Mean Squared Error (MSE): {gb_mse} ({gb_mse_percent:.2f}%)\")\n",
    "print(f\"Gradient Boosting - Root Mean Squared Error (RMSE): {gb_rmse} ({gb_rmse_percent:.2f}%)\")\n",
    "\n",
    "# Add GB predictions to the results dataframe\n",
    "results_df['Predicted_GB'] = np.nan\n",
    "results_df.loc[X_test.index, 'Predicted_GB'] = gb_predictions\n",
    "\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.plot(results_df['TIME_PERIOD'], y, label='Actual Values', color='b', linewidth = 0.8)\n",
    "plt.plot(results_df['TIME_PERIOD'], results_df['Predicted_GB'], label='GB Predicted Values', color = 'r', linestyle = '--', marker = 'o', markersize = 1)\n",
    "plt.xlabel('Time')\n",
    "plt.ylabel(dependent_var)\n",
    "plt.title(f'Actual vs Gradient Boosting Predicted Values for {dependent_var}')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "importances_gb = gb_model.feature_importances_\n",
    "indices_gb = np.argsort(importances_gb)[::-1]\n",
    "features_gb = [independent_vars[i] for i in indices_gb]\n",
    "\n",
    "# Print feature importances for Gradient Boosting\n",
    "print(\"\\nFeature Importances based on Gini Index (Gradient Boosting):\")\n",
    "for feature_name, importance in zip(features_gb, importances_gb[indices_gb]):\n",
    "    print(f\"{feature_name}: {importance:.4f}\")\n",
    "\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.bar(features_gb, importances_gb[indices_gb], color='blue')\n",
    "plt.axhline(y=0.05, color='red', linestyle='--', label='Significance Level (0.05)')\n",
    "plt.xlabel('Variables')\n",
    "plt.ylabel('Importance')\n",
    "plt.title('Variable Importance Based on Gradient Boosting')\n",
    "plt.xticks(rotation=45)\n",
    "plt.legend()  # Add legend to show the significance level\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "#-----------------------------------\n",
    "\n",
    "# Define the scoring metrics to use in cross-validation\n",
    "scoring = {'MAE': make_scorer(mean_absolute_error),\n",
    "           'MSE': make_scorer(mean_squared_error),\n",
    "           'RMSE': make_scorer(lambda y_true, y_pred: np.sqrt(mean_squared_error(y_true, y_pred)))}\n",
    "\n",
    "# Linear Regression Cross-Validation (10-fold)\n",
    "print(\"\\nLinear Regression Cross-Validation (10-fold):\")\n",
    "cv_scores_lr = cross_val_score(model_lr, X_scaled, y, cv=10, scoring='neg_mean_squared_error')\n",
    "cv_mae_lr = cross_val_score(model_lr, X_scaled, y, cv=10, scoring='neg_mean_absolute_error')\n",
    "\n",
    "# Convert negative scores to positive for interpretation\n",
    "cv_mse_lr = -cv_scores_lr\n",
    "cv_mae_lr = -cv_mae_lr\n",
    "cv_rmse_lr = np.sqrt(cv_mse_lr)\n",
    "\n",
    "# Calculate average and standard deviation for the metrics\n",
    "print(f\"Linear Regression - Average MSE: {cv_mse_lr.mean():.4f} ± {cv_mse_lr.std():.4f}\")\n",
    "print(f\"Linear Regression - Average MAE: {cv_mae_lr.mean():.4f} ± {cv_mae_lr.std():.4f}\")\n",
    "print(f\"Linear Regression - Average RMSE: {cv_rmse_lr.mean():.4f} ± {cv_rmse_lr.std():.4f}\")\n",
    "\n",
    "# Random Forest Cross-Validation (10-fold)\n",
    "print(\"\\nRandom Forest Cross-Validation (10-fold):\")\n",
    "cv_scores_rf = cross_val_score(rf_model, X_scaled, y, cv=10, scoring='neg_mean_squared_error')\n",
    "cv_mae_rf = cross_val_score(rf_model, X_scaled, y, cv=10, scoring='neg_mean_absolute_error')\n",
    "\n",
    "# Convert negative scores to positive for interpretation\n",
    "cv_mse_rf = -cv_scores_rf\n",
    "cv_mae_rf = -cv_mae_rf\n",
    "cv_rmse_rf = np.sqrt(cv_mse_rf)\n",
    "\n",
    "# Calculate average and standard deviation for the metrics\n",
    "print(f\"Random Forest - Average MSE: {cv_mse_rf.mean():.4f} ± {cv_mse_rf.std():.4f}\")\n",
    "print(f\"Random Forest - Average MAE: {cv_mae_rf.mean():.4f} ± {cv_mae_rf.std():.4f}\")\n",
    "print(f\"Random Forest - Average RMSE: {cv_rmse_rf.mean():.4f} ± {cv_rmse_rf.std():.4f}\")\n",
    "\n",
    "# Gradient Boosting Cross-Validation (10-fold)\n",
    "print(\"\\nGradient Boosting Cross-Validation (10-fold):\")\n",
    "cv_scores_gb = cross_val_score(gb_model, X_scaled, y, cv=10, scoring='neg_mean_squared_error')\n",
    "cv_mae_gb = cross_val_score(gb_model, X_scaled, y, cv=10, scoring='neg_mean_absolute_error')\n",
    "\n",
    "# Convert negative scores to positive for interpretation\n",
    "cv_mse_gb = -cv_scores_gb\n",
    "cv_mae_gb = -cv_mae_gb\n",
    "cv_rmse_gb = np.sqrt(cv_mse_gb)\n",
    "\n",
    "# Calculate average and standard deviation for the metrics\n",
    "print(f\"Gradient Boosting - Average MSE: {cv_mse_gb.mean():.4f} ± {cv_mse_gb.std():.4f}\")\n",
    "print(f\"Gradient Boosting - Average MAE: {cv_mae_gb.mean():.4f} ± {cv_mae_gb.std():.4f}\")\n",
    "print(f\"Gradient Boosting - Average RMSE: {cv_rmse_gb.mean():.4f} ± {cv_rmse_gb.std():.4f}\")\n"
   ],
   "id": "d24dab5463ca1911",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "# MACHINE LEARNING",
   "id": "1f384ea31b063031"
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "# ***Classification Analysis***",
   "id": "999203b5ed418984"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.metrics import accuracy_score, classification_report, confusion_matrix\n",
    "from sklearn.ensemble import RandomForestClassifier\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.model_selection import train_test_split\n",
    "from sklearn.model_selection import cross_val_score\n",
    "\n",
    "# Read the dataset\n",
    "dataset_file_path = f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv'\n",
    "df = pd.read_csv(dataset_file_path)\n",
    "\n",
    "# Convert the TIME_PERIOD column to datetime\n",
    "df['TIME_PERIOD'] = pd.to_datetime(df['TIME_PERIOD'])\n",
    "\n",
    "# Define dependent and independent variables for Team-1\n",
    "#dependent_var = f'UN_{country}_T_TOT'  # Replace with your dependent variable\n",
    "#independent_vars = ['Covid19_Trend_Binary', 'ECB_EXCHANGE_RATES_USD_EUR', f'Covid_{country}_Cumulative_Deaths', 'Oil_Close_Price', 'Corn_Close_Price']  # Replace with your independent variables\n",
    "\n",
    "# Define dependent and independent variables for Team-3\n",
    "dependent_var = f'INFLATION_{country}_VALUE'  # Replace with your dependent variable\n",
    "independent_vars = ['Russia_Ukraine_War_VideoCount_Binary', 'Oil_Close_Price', f'Covid_{country}_Cumulative_Deaths', 'ECB_EXCHANGE_RATES_CNY_EUR', f'BS_{country}_FS_NY', f'UN_{country}_T_TOT', 'Covid19_Trend_Binary', 'ECB_EXCHANGE_RATES_USD_EUR']  # Replace with your independent variables\n",
    "\n",
    "# Create X and y\n",
    "X = df[independent_vars]\n",
    "y = df[dependent_var]\n",
    "\n",
    "# Z-score normalization\n",
    "scaler_X = StandardScaler()\n",
    "X_scaled = scaler_X.fit_transform(X)\n",
    "\n",
    "# Check if the target variable is already binary (0 and 1)\n",
    "if set(y.unique()) == {0, 1}:\n",
    "    y_class = y\n",
    "else:\n",
    "    # Create a binary classification for the dependent variable\n",
    "    median_value = y.median()\n",
    "    y_class = (y > median_value).astype(int)  # 1 if the value is above the median, otherwise 0\n",
    "    print(f\"Class 1: > {median_value}, Class 0: ≤ {median_value}\")\n",
    "\n",
    "# Split the data into training and testing sets\n",
    "X_train, X_test, y_train, y_test = train_test_split(X_scaled, y_class, test_size=0.2, random_state=42)\n",
    "\n",
    "# Logistic Regression Model\n",
    "model_lr = LogisticRegression(random_state=42)\n",
    "model_lr.fit(X_train, y_train)\n",
    "\n",
    "# Predict the classes with Logistic Regression\n",
    "y_pred_class_lr = model_lr.predict(X_test)\n",
    "\n",
    "# Evaluate Logistic Regression\n",
    "accuracy_lr = accuracy_score(y_test, y_pred_class_lr)\n",
    "print(f'Logistic Regression Classification Accuracy: {accuracy_lr}')\n",
    "print('Logistic Regression Classification Report:')\n",
    "print(classification_report(y_test, y_pred_class_lr))\n",
    "print('Logistic Regression Confusion Matrix:')\n",
    "print(confusion_matrix(y_test, y_pred_class_lr))\n",
    "\n",
    "# Visualize the results of Logistic Regression\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.scatter(df.iloc[y_test.index]['TIME_PERIOD'], y_test, label='Actual Class', color='b', alpha=0.6)\n",
    "plt.scatter(df.iloc[y_test.index]['TIME_PERIOD'], y_pred_class_lr, label='Predicted Class', color='r', alpha=0.6, marker='x')\n",
    "plt.xlabel('Time')\n",
    "plt.ylabel('Class')\n",
    "if 'median_value' in locals():\n",
    "    plt.title(f'Actual vs Predicted Classes for {dependent_var} - Logistic Regression\\nClass 1: > {median_value}, Class 0: ≤ {median_value}')\n",
    "else:\n",
    "    plt.title(f'Actual vs Predicted Classes for {dependent_var} - Logistic Regression')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "# Feature importance based on Logistic Regression coefficients\n",
    "coefficients = model_lr.coef_[0]\n",
    "indices = np.argsort(abs(coefficients))[::-1]\n",
    "\n",
    "# Print the feature ranking based on the absolute value of coefficients\n",
    "print(\"Feature ranking based on Logistic Regression coefficients (absolute value):\")\n",
    "\n",
    "for f in range(X.shape[1]):\n",
    "    print(f\"{f + 1}. feature {independent_vars[indices[f]]} (coefficient: {coefficients[indices[f]]})\")\n",
    "\n",
    "# Plot the feature importances based on Logistic Regression coefficients\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.title(\"Feature Importance Based on Logistic Regression Coefficients\")\n",
    "plt.bar(range(X.shape[1]), coefficients[indices], align=\"center\")\n",
    "plt.xticks(range(X.shape[1]), [independent_vars[i] for i in indices], rotation=90)\n",
    "plt.xlim([-1, X.shape[1]])\n",
    "plt.ylabel('Coefficient Value')\n",
    "plt.tight_layout()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "# Random Forest Model\n",
    "model_rf = RandomForestClassifier(n_estimators=100, random_state=42)\n",
    "model_rf.fit(X_train, y_train)\n",
    "\n",
    "# Predict the classes with Random Forest\n",
    "y_pred_class_rf = model_rf.predict(X_test)\n",
    "\n",
    "# Evaluate Random Forest\n",
    "accuracy_rf = accuracy_score(y_test, y_pred_class_rf)\n",
    "print(f'Random Forest Classification Accuracy: {accuracy_rf}')\n",
    "print('Random Forest Classification Report:')\n",
    "print(classification_report(y_test, y_pred_class_rf))\n",
    "print('Random Forest Confusion Matrix:')\n",
    "print(confusion_matrix(y_test, y_pred_class_rf))\n",
    "\n",
    "# Visualize the results of Random Forest\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.scatter(df.iloc[y_test.index]['TIME_PERIOD'], y_test, label='Actual Class', color='b', alpha=0.6)\n",
    "plt.scatter(df.iloc[y_test.index]['TIME_PERIOD'], y_pred_class_rf, label='Predicted Class', color='r', alpha=0.6, marker='x')\n",
    "plt.xlabel('Time')\n",
    "plt.ylabel('Class')\n",
    "if 'median_value' in locals():\n",
    "    plt.title(f'Actual vs Predicted Classes for {dependent_var} - Random Forest Classification\\nClass 1: > {median_value}, Class 0: ≤ {median_value}')\n",
    "else:\n",
    "    plt.title(f'Actual vs Predicted Classes for {dependent_var} - Random Forest Classification')\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "# Feature importance based on Gini index from Random Forest\n",
    "importances = model_rf.feature_importances_\n",
    "indices = np.argsort(importances)[::-1]\n",
    "\n",
    "# Print the feature ranking\n",
    "print(\"Feature ranking based on Gini index:\")\n",
    "\n",
    "for f in range(X.shape[1]):\n",
    "    print(f\"{f + 1}. feature {independent_vars[indices[f]]} ({importances[indices[f]]})\")\n",
    "\n",
    "# Plot the feature importances of the forest\n",
    "plt.figure(figsize=(10, 6))\n",
    "plt.title(\"Feature Importance Based on Gini Index (Random Forest)\")\n",
    "plt.bar(range(X.shape[1]), importances[indices], align=\"center\")\n",
    "plt.xticks(range(X.shape[1]), [independent_vars[i] for i in indices], rotation=90)\n",
    "plt.xlim([-1, X.shape[1]])\n",
    "plt.ylabel('Importance')\n",
    "plt.tight_layout()\n",
    "plt.grid(True)\n",
    "plt.show()\n",
    "\n",
    "#-----------------------------\n",
    "\n",
    "# Logistic Regression Model with Cross-Validation\n",
    "model_lr_cv = LogisticRegression(random_state=42)\n",
    "cv_scores_lr = cross_val_score(model_lr_cv, X_scaled, y_class, cv=10)  # 10-fold cross-validation\n",
    "\n",
    "# Print Logistic Regression Cross-Validation Results\n",
    "print(f'Logistic Regression Cross-Validation Accuracy: {cv_scores_lr.mean()} ± {cv_scores_lr.std()}')\n",
    "\n",
    "# Random Forest Model with Cross-Validation\n",
    "model_rf_cv = RandomForestClassifier(n_estimators=100, random_state=42)\n",
    "cv_scores_rf = cross_val_score(model_rf_cv, X_scaled, y_class, cv=10)  # 10-fold cross-validation\n",
    "\n",
    "# Print Random Forest Cross-Validation Results\n",
    "print(f'Random Forest Cross-Validation Accuracy: {cv_scores_rf.mean()} ± {cv_scores_rf.std()}')\n",
    "\n",
    "# Train Logistic Regression on the entire dataset for feature importance\n",
    "model_lr_cv.fit(X_scaled, y_class)\n",
    "\n",
    "# Feature importance based on Logistic Regression coefficients\n",
    "coefficients_cv = model_lr_cv.coef_[0]\n",
    "indices_lr_cv = np.argsort(abs(coefficients_cv))[::-1]\n",
    "\n",
    "# Print the feature ranking based on Logistic Regression coefficients with cross-validation\n",
    "print(\"Feature ranking based on Logistic Regression coefficients (absolute value) with cross-validation:\")\n",
    "for f in range(X.shape[1]):\n",
    "    print(f\"{f + 1}. feature {independent_vars[indices_lr_cv[f]]} (coefficient: {coefficients_cv[indices_lr_cv[f]]})\")\n",
    "\n",
    "# Train Random Forest on the entire dataset for feature importance\n",
    "model_rf_cv.fit(X_scaled, y_class)\n",
    "\n",
    "# Feature importance based on Gini index from Random Forest with cross-validation\n",
    "importances_rf_cv = model_rf_cv.feature_importances_\n",
    "indices_rf_cv = np.argsort(importances_rf_cv)[::-1]\n",
    "\n",
    "# Print the feature ranking based on Gini index with cross-validation\n",
    "print(\"Feature ranking based on Gini index (Random Forest) with cross-validation:\")\n",
    "for f in range(X.shape[1]):\n",
    "    print(f\"{f + 1}. feature {independent_vars[indices_rf_cv[f]]} ({importances_rf_cv[indices_rf_cv[f]]})\")\n"
   ],
   "id": "3603d12cb3f80baa",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "# ***Clustering Analysis***",
   "id": "3df5220af3eb9232"
  },
  {
   "metadata": {},
   "cell_type": "code",
   "source": [
    "import pandas as pd\n",
    "from sklearn.preprocessing import StandardScaler, MinMaxScaler\n",
    "from sklearn.cluster import KMeans, AgglomerativeClustering\n",
    "import matplotlib.pyplot as plt\n",
    "from sklearn.metrics import silhouette_score\n",
    "from scipy.cluster.hierarchy import dendrogram, linkage\n",
    "\n",
    "# Read the dataset\n",
    "dataset_file_path = f'EDA_ML_DATA/{country}/updated_merged_data_{country}.csv'\n",
    "df = pd.read_csv(dataset_file_path)\n",
    "\n",
    "# Define the new sets of features\n",
    "economic_unemployment_features = [\n",
    "    'ESTR_EU000A2QQF08_CI', f'CISS_{country}_SS_CIN', f'GOVDEBT_{country}', #'ESTR_EU000A2X2A25_TT',\n",
    "    f'INFLATION_{country}_VALUE', f'GDP_GROWTH_{country}',  \n",
    "    f'UN_{country}_T_TOT'\n",
    "]\n",
    "\n",
    "market_risk_and_confidence_features = [\n",
    "    f'BS_{country}_CSMCI', f'BS_{country}_FS_NY', f'BS_{country}_GES_NY', f'BS_{country}_MP_NY', \n",
    "    f'BS_{country}_PT_NY', f'BS_{country}_SV_NY', 'ECB_EXCHANGE_RATES_USD_EUR',\n",
    "    'ECB_EXCHANGE_RATES_CNY_EUR', 'Oil_Close_Price', 'Gas_Close_Price'    \n",
    "]\n",
    "\n",
    "# Function to perform clustering and visualize results\n",
    "def perform_clustering(features, feature_set_name, n_clusters):\n",
    "    # Create a new DataFrame with only the selected features\n",
    "    df_features = df[features]\n",
    "\n",
    "    # Normalize the data with StandardScaler\n",
    "    scaler = StandardScaler()\n",
    "    df_features_scaled = scaler.fit_transform(df_features)\n",
    "    \n",
    "    # Elbow Method and Silhouette Scores for determining the optimal number of clusters\n",
    "    inertia = []\n",
    "    silhouette_scores = []\n",
    "    K = range(2, 11)\n",
    "\n",
    "    for k in K:\n",
    "        kmeans = KMeans(n_clusters=k, random_state=42)\n",
    "        kmeans.fit(df_features_scaled)\n",
    "        inertia.append(kmeans.inertia_)\n",
    "        silhouette_scores.append(silhouette_score(df_features_scaled, kmeans.labels_))\n",
    "\n",
    "    # Plot the Elbow Method\n",
    "    plt.figure(figsize=(10, 6))\n",
    "    plt.plot(K, inertia, 'bo-')\n",
    "    plt.xlabel('Number of clusters (k)')\n",
    "    plt.ylabel('Inertia')\n",
    "    plt.title(f'Elbow Method for Optimal Number of Clusters ({feature_set_name})')\n",
    "    plt.grid(True)\n",
    "    plt.show()\n",
    "\n",
    "    # Plot the Silhouette Scores\n",
    "    plt.figure(figsize=(10, 6))\n",
    "    plt.plot(K, silhouette_scores, 'bo-')\n",
    "    plt.xlabel('Number of clusters (k)')\n",
    "    plt.ylabel('Silhouette Score')\n",
    "    plt.title(f'Silhouette Scores for Optimal Number of Clusters ({feature_set_name})')\n",
    "    plt.grid(True)\n",
    "    plt.show()\n",
    "\n",
    "    # Apply K-means clustering with the specified number of clusters\n",
    "    kmeans = KMeans(n_clusters=n_clusters, random_state=42)\n",
    "    df[f'KMeans_Cluster_{feature_set_name}'] = kmeans.fit_predict(df_features_scaled)\n",
    "\n",
    "    # Visualize the K-means clustering results using PCA\n",
    "    from sklearn.decomposition import PCA\n",
    "\n",
    "    pca = PCA(n_components=2)\n",
    "    df_pca = pca.fit_transform(df_features_scaled)\n",
    "\n",
    "    plt.figure(figsize=(10, 6))\n",
    "    plt.scatter(df_pca[:, 0], df_pca[:, 1], c=df[f'KMeans_Cluster_{feature_set_name}'], cmap='viridis')\n",
    "    plt.xlabel('Principal Component 1')\n",
    "    plt.ylabel('Principal Component 2')\n",
    "    plt.title(f'K-means Clustering of {feature_set_name}')\n",
    "    plt.colorbar(label='Cluster')\n",
    "    plt.grid(True)\n",
    "    plt.show()\n",
    "\n",
    "    # Analyze clusters by calculating mean values of the features within each cluster\n",
    "    kmeans_cluster_means = df.groupby(f'KMeans_Cluster_{feature_set_name}')[features].mean()\n",
    "\n",
    "    # Print mean values for interpretation\n",
    "    print(f\"\\nMean values of {feature_set_name} for each K-means cluster:\")\n",
    "    print(kmeans_cluster_means)\n",
    "\n",
    "    # Visualize the mean values to help with interpretation\n",
    "    #fig, ax = plt.subplots(figsize=(16, 10))  # Adjusting the figsize to make the plot larger\n",
    "    #kmeans_cluster_means.plot(kind='bar', ax=ax, title=f'Mean Values of {feature_set_name} by K-means Cluster', logy=True)\n",
    "    #ax.legend(bbox_to_anchor=(1.05, 1), loc='upper left')  # Move the legend outside of the plot\n",
    "    #plt.tight_layout()\n",
    "    #plt.show()\n",
    "    \n",
    "    # Visualize the mean values to help with interpretation\n",
    "    fig, ax = plt.subplots(figsize=(16, 10))  # Adjusting the figsize to make the plot larger\n",
    "    kmeans_cluster_means.plot(kind='bar', ax=ax, title=f'Mean Values of {feature_set_name} by K-means Cluster')\n",
    "    \n",
    "    # Προσαρμογή του εύρους του άξονα y για να περιλαμβάνει αρνητικές και θετικές τιμές\n",
    "    ax.set_ylim([kmeans_cluster_means.values.min() - 1, kmeans_cluster_means.values.max() + 1])\n",
    "    \n",
    "    ax.legend(bbox_to_anchor=(1.05, 1), loc='upper left')  # Move the legend outside of the plot\n",
    "    plt.tight_layout()\n",
    "    plt.show()\n",
    "\n",
    "\n",
    "# Define the number of clusters for each feature set\n",
    "economic_unemployment_clusters = 4  # Change this value to the desired number of clusters for economic and unemployment features\n",
    "market_risk_confidence_clusters = 3  # Change this value to the desired number of clusters for market risk and confidence features\n",
    "\n",
    "# Perform clustering on economic and unemployment features\n",
    "perform_clustering(economic_unemployment_features, 'Economic and Unemployment Features', economic_unemployment_clusters)\n",
    "\n",
    "# Perform clustering on market risk features\n",
    "perform_clustering(market_risk_and_confidence_features, 'Market Risk and Confidence Features', market_risk_confidence_clusters)\n",
    "\n",
    "\n",
    "#---------------------------\n",
    "\n",
    "# Function to perform Agglomerative Clustering and visualize results\n",
    "def perform_agglomerative_clustering(features, feature_set_name, n_clusters):\n",
    "    # Create a new DataFrame with only the selected features\n",
    "    df_features = df[features]\n",
    "\n",
    "    # Normalize the data with StandardScaler\n",
    "    scaler = StandardScaler()\n",
    "    df_features_scaled = scaler.fit_transform(df_features)\n",
    "    \n",
    "    # Normalize the data with Min-Max Scaler\n",
    "    #scaler = MinMaxScaler()\n",
    "    #df_features_scaled = scaler.fit_transform(df_features)\n",
    "\n",
    "    # Apply Agglomerative Clustering\n",
    "    agglomerative = AgglomerativeClustering(n_clusters=n_clusters)\n",
    "    df[f'Agglomerative_Cluster_{feature_set_name}'] = agglomerative.fit_predict(df_features_scaled)\n",
    "\n",
    "    # Visualize the Agglomerative Clustering results using PCA\n",
    "    from sklearn.decomposition import PCA\n",
    "\n",
    "    pca = PCA(n_components=2)\n",
    "    df_pca = pca.fit_transform(df_features_scaled)\n",
    "\n",
    "    plt.figure(figsize=(10, 6))\n",
    "    plt.scatter(df_pca[:, 0], df_pca[:, 1], c=df[f'Agglomerative_Cluster_{feature_set_name}'], cmap='viridis')\n",
    "    plt.xlabel('Principal Component 1')\n",
    "    plt.ylabel('Principal Component 2')\n",
    "    plt.title(f'Agglomerative Clustering of {feature_set_name}')\n",
    "    plt.colorbar(label='Cluster')\n",
    "    plt.grid(True)\n",
    "    plt.show()\n",
    "\n",
    "    # Generate and plot the dendrogram\n",
    "    Z = linkage(df_features_scaled, 'ward')\n",
    "    plt.figure(figsize=(10, 6))\n",
    "    dendrogram(Z)\n",
    "    plt.title(f'Dendrogram for {feature_set_name}')\n",
    "    plt.xlabel('Samples')\n",
    "    plt.ylabel('Distance')\n",
    "    plt.show()\n",
    "\n",
    "    # Analyze clusters by calculating mean values of the features within each cluster\n",
    "    agglomerative_cluster_means = df.groupby(f'Agglomerative_Cluster_{feature_set_name}')[features].mean()\n",
    "\n",
    "    # Print mean values for interpretation\n",
    "    print(f\"\\nMean values of {feature_set_name} for each Agglomerative cluster:\")\n",
    "    print(agglomerative_cluster_means)\n",
    "\n",
    "    # Visualize the mean values to help with interpretation\n",
    "    #fig, ax = plt.subplots(figsize=(16, 10))  # Adjusting the figsize to make the plot larger\n",
    "    #agglomerative_cluster_means.plot(kind='bar', ax=ax, title=f'Mean Values of {feature_set_name} by Agglomerative Cluster', logy=True)\n",
    "    #ax.legend(bbox_to_anchor=(1.05, 1), loc='upper left')  # Move the legend outside of the plot\n",
    "    #plt.tight_layout()\n",
    "    #plt.show()\n",
    "    \n",
    "    # Visualize the mean values to help with interpretation\n",
    "    fig, ax = plt.subplots(figsize=(16, 10))  # Adjusting the figsize to make the plot larger\n",
    "    agglomerative_cluster_means.plot(kind='bar', ax=ax, title=f'Mean Values of {feature_set_name} by Agglomerative Cluster')\n",
    "    \n",
    "    ax.set_ylim([agglomerative_cluster_means.values.min() - 1, agglomerative_cluster_means.values.max() + 1])\n",
    "    \n",
    "    ax.legend(bbox_to_anchor=(1.05, 1), loc='upper left')  # Move the legend outside of the plot\n",
    "    plt.tight_layout()\n",
    "    plt.show()\n",
    "    \n",
    "    # Define the number of clusters for each feature set\n",
    "economic_unemployment_clusters = 4  # Change this value to the desired number of clusters for economic and unemployment features\n",
    "market_risk_confidence_clusters = 3  # Change this value to the desired number of clusters for market risk and confidence features\n",
    "\n",
    "# Perform Agglomerative Clustering on economic and unemployment features\n",
    "perform_agglomerative_clustering(economic_unemployment_features, 'Economic and Unemployment Features', economic_unemployment_clusters)\n",
    "\n",
    "# Perform Agglomerative Clustering on market risk features\n",
    "perform_agglomerative_clustering(market_risk_and_confidence_features, 'Market Risk and Confidence Features', market_risk_confidence_clusters)\n"
   ],
   "id": "a11d70c069346376",
   "outputs": [],
   "execution_count": null
  },
  {
   "metadata": {},
   "cell_type": "code",
   "outputs": [],
   "execution_count": null,
   "source": "",
   "id": "7a7c52a9253978b8"
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}