RARE-JupyterNotebook.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "529e879c",
   "metadata": {},
   "source": [
    "# RARE\n",
    "RARE: Relevant Association Rare-variant-bin Evolver is an evolutionary algorithm approach to binning rare variants as a rare variant association analysis tool\n",
    "## What is RARE? \n",
    "RARE is an evolutionary algorithm that constructs bins of rare variant features with relevant association to class (univariate and/or multivariate interactions) through the following steps:\n",
    "1. Random bin initializaiton or expert knowledge input\n",
    "2. Repeated evolutionary cycles consisting of:\n",
    "    - Candidate bin evaluation with univariate scoring (chi-square test) or Relief-based scoring (MultiSURF algorithm); note: new scoring options currently under testing\n",
    "    - Genetic operations (parent selection, crossover, and mutation) to generate the next generation of candidate bins\n",
    "3. Final bin evaluation and summary of top bins\n",
    "\n",
    "### What this notebook contains:  \n",
    "1. An example of a dataset RARE can be run on\n",
    "2. How to run RARE \n",
    "3. How to interpret RARE's results"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "e0fbd137",
   "metadata": {},
   "outputs": [],
   "source": [
    "from skrare.methods import *\n",
    "from IPython.display import display\n",
    "import warnings\n",
    "warnings.filterwarnings('ignore')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "cb4b9edc",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "sim_data = pd.read_csv('Experiment1.csv')"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2ad2f7a7",
   "metadata": {},
   "source": [
    "For this example, the data was created with 1000 instances, 50 total rare variant features, and 10 predictive rare variant features. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "8f3325b7",
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div>\n",
       "<style scoped>\n",
       "    .dataframe tbody tr th:only-of-type {\n",
       "        vertical-align: middle;\n",
       "    }\n",
       "\n",
       "    .dataframe tbody tr th {\n",
       "        vertical-align: top;\n",
       "    }\n",
       "\n",
       "    .dataframe thead th {\n",
       "        text-align: right;\n",
       "    }\n",
       "</style>\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>P_1</th>\n",
       "      <th>P_2</th>\n",
       "      <th>P_3</th>\n",
       "      <th>P_4</th>\n",
       "      <th>P_5</th>\n",
       "      <th>P_6</th>\n",
       "      <th>P_7</th>\n",
       "      <th>P_8</th>\n",
       "      <th>P_9</th>\n",
       "      <th>P_10</th>\n",
       "      <th>...</th>\n",
       "      <th>R_32</th>\n",
       "      <th>R_33</th>\n",
       "      <th>R_34</th>\n",
       "      <th>R_35</th>\n",
       "      <th>R_36</th>\n",
       "      <th>R_37</th>\n",
       "      <th>R_38</th>\n",
       "      <th>R_39</th>\n",
       "      <th>R_40</th>\n",
       "      <th>Class</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>...</th>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "      <td>...</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>995</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>996</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>997</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>2.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>998</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>999</th>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>...</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>0.0</td>\n",
       "      <td>1.0</td>\n",
       "      <td>0.0</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "<p>1000 rows × 51 columns</p>\n",
       "</div>"
      ],
      "text/plain": [
       "     P_1  P_2  P_3  P_4  P_5  P_6  P_7  P_8  P_9  P_10  ...  R_32  R_33  R_34  \\\n",
       "0    0.0  0.0  0.0  0.0  0.0  0.0  0.0  1.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "1    0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "2    0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "3    0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "4    0.0  0.0  0.0  1.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   1.0   \n",
       "..   ...  ...  ...  ...  ...  ...  ...  ...  ...   ...  ...   ...   ...   ...   \n",
       "995  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "996  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "997  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "998  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "999  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0  0.0   0.0  ...   0.0   0.0   0.0   \n",
       "\n",
       "     R_35  R_36  R_37  R_38  R_39  R_40  Class  \n",
       "0     0.0   0.0   0.0   0.0   0.0   1.0    1.0  \n",
       "1     0.0   0.0   0.0   0.0   0.0   0.0    0.0  \n",
       "2     0.0   0.0   0.0   0.0   1.0   0.0    0.0  \n",
       "3     0.0   0.0   0.0   0.0   0.0   0.0    0.0  \n",
       "4     0.0   0.0   0.0   0.0   0.0   0.0    1.0  \n",
       "..    ...   ...   ...   ...   ...   ...    ...  \n",
       "995   0.0   0.0   0.0   0.0   0.0   0.0    0.0  \n",
       "996   0.0   0.0   0.0   1.0   0.0   0.0    0.0  \n",
       "997   0.0   2.0   0.0   0.0   0.0   0.0    0.0  \n",
       "998   0.0   0.0   0.0   0.0   0.0   0.0    0.0  \n",
       "999   0.0   0.0   0.0   0.0   0.0   1.0    0.0  \n",
       "\n",
       "[1000 rows x 51 columns]"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "display(sim_data)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "48b9b0e6",
   "metadata": {},
   "source": [
    "### Parameters for RARE:\n",
    "1. given_starting_point: whether or not expert knowledge is being inputted (True or False)\n",
    "2. amino_acid_start_point: if RARE is starting with expert knowledge, input the list of features here; otherwise None\n",
    "3. amino_acid_bins_start_point: if RARE is starting with expert knowledge, input the list of bins of features here; otherwise None\n",
    "4. iterations: the number of evolutionary cycles RARE will run\n",
    "5. original_feature_matrix: the dataset\n",
    "6. label_name: label for the class/endpoint column in the dataset (e.g., 'Class')\n",
    "7. rare_variant_MAF_cutoff: the minor allele frequency cutoff separating common features from rare variant features\n",
    "8. set_number_of_bins: the population size of candidate bins\n",
    "9. min_features_per_group: the minimum number of features in a bin\n",
    "10. max_number_of_groups_with_feature: the maximum number of bins containing a feature\n",
    "11. scoring_method: 'Univariate', 'Relief', or 'Relief only on bin and common features'\n",
    "12. score_based_on_sample: if Relief scoring is used, whether or not bin evaluation is done based on a sample of instances rather than the whole dataset\n",
    "13. score_with_common_variables: if Relief scoring is used, whether or not common features should be used as context for evaluating rare variant bins\n",
    "14. instance_sample_size: if bin evaluation is done based on a sample of instances, input the sample size here\n",
    "15. crossover_probability: the probability of each feature in an offspring bin to crossover to the paired offspring bin (recommendation: 0.5 to 0.8)\n",
    "16. mutation_probability: the probability of each feature in a bin to be deleted (a proportionate probability is automatically applied on each feature outside the bin to be added (recommendation: 0.05 to 0.5 depending on situation and number of iterations run)\n",
    "17. elitism_parameter: the proportion of elite bins in the current generation to be preserved for the next evolutionary cycle (recommendation: 0.2 to 0.8 depending on conservativeness of approach and number of iterations run)\n",
    "18. random_seed: the seed value needed to generate a random number to create randomness in crossover and mutation\n",
    "19. bin_size_variability_constraint: sets the max bin size of children to be n times the size of their sibling (recommendation: 2, with larger or smaller values the population would trend heavily towards small or large bins without exploring the search space)\n",
    "20. max_features_per_bin: sets a max value for the number of features per bin"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f9424cba",
   "metadata": {
    "scrolled": false
   },
   "outputs": [],
   "source": [
    "%%capture\n",
    "for replicate in range (0,1):\n",
    "    print('Experiment:');\n",
    "    bin_feature_matrix, common_features_and_bins_matrix, amino_acid_bins, amino_acid_bin_scores, rare_feature_MAF_dict, common_feature_MAF_dict, rare_feature_df, common_feature_df, MAF_0_features = RARE_v2(False, None, None, 100, sim_data, 'Class', 0.05, 50, 5, 25, 'Relief', True, False, 500, 0.8, 0.1, 0.4, 10, 2, None);\n",
    "        \n",
    "    # Summary of the best bins\n",
    "    summary = Top_Rare_Variant_Bins_Summary(rare_feature_df, 'Class', amino_acid_bins, amino_acid_bin_scores, rare_feature_MAF_dict, 5, bin_feature_matrix)\n",
    "    \n",
    "    #Exporting summary to a CSV file\n",
    "    summary.to_csv('Experiment1Results.csv', index=False)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2f097045",
   "metadata": {},
   "source": [
    "### Results\n",
    "The bins are ranked in order of their MultiSURF or Univariate Score with higher scores corresponding to better bins. In this example, we're performing relief-based scoring using the MultiSURF Algorithm. This algorithm allows us to evaluate a given bin relative to all other bins in the population, selecting features in a way that is sensitive to feature interactions. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "50e13a1f",
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "pd.set_option('display.max_rows', None)\n",
    "pd.set_option('display.max_columns', None)\n",
    "pd.set_option('display.width', None)\n",
    "pd.set_option('display.max_colwidth', None)\n",
    "display(summary)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "4225bb9a",
   "metadata": {},
   "source": [
    "# SciKit Learn Package\n",
    "RARE has been implemented as a Scikit-learn package. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "abadd572",
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv('Experiment1.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "outputs": [],
   "source": [
    "from sklearn.pipeline import Pipeline\n",
    "from skrare.rare import RARE"
   ],
   "metadata": {
    "collapsed": false
   }
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "0bc7b03d",
   "metadata": {},
   "outputs": [],
   "source": [
    "pipe_RARE = Pipeline(steps=[(\"RARE\", RARE(given_starting_point=False, amino_acid_start_point=None, \n",
    "amino_acid_bins_start_point=None, iterations=1000, label_name= “Class”, rare_variant_MAF_cutoff=0.05, set_number_of_bins=50, min_features_per_group=5, max_number_of_groups_with_feature=25, scoring_method=‘Relief’, score_based_on_sample=True, score_with_common_variables=False, instance_sample_size=500, crossover_probability=0.8, mutation_probability=0.1, elitism_parameter=0.4, random_seed=None, bin_size_variability=None))])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "7fe47163",
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}