{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"#import libraiares:\nimport tensorflow as tf\nimport pandas as pd\nimport numpy as np\n\n# we will use bio python to read the sequences\nfrom Bio import SeqIO\n\nimport seaborn as sns\nimport matplotlib.pyplot as plt\n\n\nfrom sklearn.impute import KNNImputer, SimpleImputer\n#import statsmodels.api as sm\n\nfrom sklearn.model_selection import train_test_split, RandomizedSearchCV, GridSearchCV\n#from sklearn.ensemble import RandomForestRegressor\nfrom scipy.stats import spearmanr\nfrom sklearn.metrics import mean_absolute_error\nfrom xgboost import XGBRegressor\nfrom sklearn.preprocessing import StandardScaler\nimport re\nfrom scipy.stats import randint\n\nimport plotly.express as px\nfrom transformers import BertModel, BertTokenizer\nimport torch\n\n# Required for progressbar widget\nimport progressbar\nimport time\nt0start = time.time() \nfrom tqdm import tqdm","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:20.691927Z","iopub.execute_input":"2023-07-12T08:04:20.692273Z","iopub.status.idle":"2023-07-12T08:04:35.866933Z","shell.execute_reply.started":"2023-07-12T08:04:20.692247Z","shell.execute_reply":"2023-07-12T08:04:35.866268Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Problem definition:\nThe goal of this competition is to predict the function of a set of proteins.\nThe accurate assignment of biological function to the protein is key to understanding life at the molecular level. However, assigning function to any specific protein can be made difficult due to the multiple functions many proteins have, along with their ability to interact with multiple partners.\n\nWe are provided with a set of protein sequences on which the participants are asked to predict Gene Ontology (GO) terms in each of the three subontologies: Molecular Function (MF), Biological Process (BP), and Cellular Component (CC). This set of sequences is referred to as test superset.\n\n# Data:\n1. train_sequences.fasta - amino acid sequences for proteins in training set\n2. testsuperset.fasta - amino acid sequences for proteins on which the predictions should be made\n3. train_terms.tsv - the training set of proteins and corresponding annotated GO terms\n4. train_taxonomy.tsv - taxon ID for proteins in training set\n5. go-basic.obo - ontology graph structure\n6. testsuperset-taxon-list.tsv - taxon ID for proteins in test superset (Note: you may need to use encoding=\"ISO-8859-1\" to read this file in pandas)\n7. IA.txt - Information Accretion for each term. This is used to weight precision and recall (see Evaluation)\n8. sample_submission.csv - a sample submission file in the correct format","metadata":{"execution":{"iopub.status.busy":"2023-06-29T08:12:34.250016Z","iopub.execute_input":"2023-06-29T08:12:34.250375Z","iopub.status.idle":"2023-06-29T08:12:34.261443Z","shell.execute_reply.started":"2023-06-29T08:12:34.250346Z","shell.execute_reply":"2023-06-29T08:12:34.260309Z"}}},{"cell_type":"code","source":"# path to the train and test fasta files\ntrain_fasta = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\ntest_fasta = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:35.868582Z","iopub.execute_input":"2023-07-12T08:04:35.869125Z","iopub.status.idle":"2023-07-12T08:04:35.874751Z","shell.execute_reply.started":"2023-07-12T08:04:35.86909Z","shell.execute_reply":"2023-07-12T08:04:35.873819Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# read train fasta files\ntrain_sequences = SeqIO.parse(train_fasta, 'fasta')","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:35.876492Z","iopub.execute_input":"2023-07-12T08:04:35.876939Z","iopub.status.idle":"2023-07-12T08:04:35.915451Z","shell.execute_reply.started":"2023-07-12T08:04:35.876906Z","shell.execute_reply":"2023-07-12T08:04:35.914605Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"sequences_dict = SeqIO.to_dict(train_sequences)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:35.917353Z","iopub.execute_input":"2023-07-12T08:04:35.917748Z","iopub.status.idle":"2023-07-12T08:04:39.095372Z","shell.execute_reply.started":"2023-07-12T08:04:35.917716Z","shell.execute_reply":"2023-07-12T08:04:39.094134Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# take a look at the first three instances:\nfirst_three_keys = list(sequences_dict.keys())[:3]\nfirst_three_instances = {key: sequences_dict[key] for key in first_three_keys}\n\n# Printing the first three instances\nfor seq_id, seq_record in first_three_instances.items():\n    print(\"ID:\", seq_id)\n    print(\"Sequence:\", seq_record.seq)\n    print('Description:', seq_record.description)\n    print()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:39.09656Z","iopub.execute_input":"2023-07-12T08:04:39.096849Z","iopub.status.idle":"2023-07-12T08:04:39.108154Z","shell.execute_reply.started":"2023-07-12T08:04:39.096826Z","shell.execute_reply":"2023-07-12T08:04:39.107566Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#create a dataframe of the above fasta files dictionary:\ndata = []\nfor seq_id, seq_record in sequences_dict.items():\n    sequence = str(seq_record.seq)\n    seq_des = str(seq_record.description)\n    data.append({\"ID\": seq_id, \"Sequence\": sequence, \"description\": seq_des})","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:39.108999Z","iopub.execute_input":"2023-07-12T08:04:39.109597Z","iopub.status.idle":"2023-07-12T08:04:39.352914Z","shell.execute_reply.started":"2023-07-12T08:04:39.109576Z","shell.execute_reply":"2023-07-12T08:04:39.351389Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_sequences = pd.DataFrame(data)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:39.3546Z","iopub.execute_input":"2023-07-12T08:04:39.354915Z","iopub.status.idle":"2023-07-12T08:04:39.499789Z","shell.execute_reply.started":"2023-07-12T08:04:39.354888Z","shell.execute_reply":"2023-07-12T08:04:39.498757Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_sequences.head()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:40.888548Z","iopub.execute_input":"2023-07-12T08:04:40.888891Z","iopub.status.idle":"2023-07-12T08:04:40.916335Z","shell.execute_reply.started":"2023-07-12T08:04:40.888865Z","shell.execute_reply":"2023-07-12T08:04:40.915329Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#create a column with length of sequences:\ndf_sequences['sequence_len'] = df_sequences['Sequence'].apply(len)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:41.96091Z","iopub.execute_input":"2023-07-12T08:04:41.961264Z","iopub.status.idle":"2023-07-12T08:04:42.01409Z","shell.execute_reply.started":"2023-07-12T08:04:41.961236Z","shell.execute_reply":"2023-07-12T08:04:42.013165Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# read test fasta files\ntest_sequences = SeqIO.parse(test_fasta, 'fasta')","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:42.67696Z","iopub.execute_input":"2023-07-12T08:04:42.6773Z","iopub.status.idle":"2023-07-12T08:04:42.682856Z","shell.execute_reply.started":"2023-07-12T08:04:42.677274Z","shell.execute_reply":"2023-07-12T08:04:42.681639Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Remove duplicates based on ID\ntest_sequences = {rec.id: rec for rec in test_sequences}.values()\n\n# Convert the sequences to a dictionary\nsequences_dict_test = SeqIO.to_dict(test_sequences)\n","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:43.288906Z","iopub.execute_input":"2023-07-12T08:04:43.289272Z","iopub.status.idle":"2023-07-12T08:04:46.109994Z","shell.execute_reply.started":"2023-07-12T08:04:43.289246Z","shell.execute_reply":"2023-07-12T08:04:46.109009Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#create a dataframe of the above fasta files dictionary:\ndata_test = []\nfor seq_id, seq_record in sequences_dict_test.items():\n    sequence = str(seq_record.seq)\n    seq_des = str(seq_record.description)\n    data_test.append({\"ID\": seq_id, \"Sequence\": sequence, \"description\": seq_des})\n\ndf_sequences_test = pd.DataFrame(data_test)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:46.112092Z","iopub.execute_input":"2023-07-12T08:04:46.112917Z","iopub.status.idle":"2023-07-12T08:04:46.447425Z","shell.execute_reply.started":"2023-07-12T08:04:46.112885Z","shell.execute_reply":"2023-07-12T08:04:46.446513Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#create a column with length of sequences:\ndf_sequences_test['sequence_len'] = df_sequences_test['Sequence'].apply(len)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:46.450592Z","iopub.execute_input":"2023-07-12T08:04:46.450875Z","iopub.status.idle":"2023-07-12T08:04:46.504474Z","shell.execute_reply.started":"2023-07-12T08:04:46.450852Z","shell.execute_reply":"2023-07-12T08:04:46.503089Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_sequences_test.head()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:47.992839Z","iopub.execute_input":"2023-07-12T08:04:47.9932Z","iopub.status.idle":"2023-07-12T08:04:48.003623Z","shell.execute_reply.started":"2023-07-12T08:04:47.993174Z","shell.execute_reply":"2023-07-12T08:04:48.002524Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Load the terms /label data:\ntrain_terms = pd.read_csv(\"/kaggle/input/cafa-5-protein-function-prediction/Train/train_terms.tsv\",sep=\"\\t\")\nprint(train_terms.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:49.864646Z","iopub.execute_input":"2023-07-12T08:04:49.864955Z","iopub.status.idle":"2023-07-12T08:04:52.604222Z","shell.execute_reply.started":"2023-07-12T08:04:49.864931Z","shell.execute_reply":"2023-07-12T08:04:52.60334Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#Load the taxanomy Dataset\ntrain_tax = pd.read_csv(\"/kaggle/input/cafa-5-protein-function-prediction/Train/train_taxonomy.tsv\",sep=\"\\t\")\nprint(train_tax.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:04:52.606565Z","iopub.execute_input":"2023-07-12T08:04:52.606919Z","iopub.status.idle":"2023-07-12T08:04:52.693838Z","shell.execute_reply.started":"2023-07-12T08:04:52.606893Z","shell.execute_reply":"2023-07-12T08:04:52.692674Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":" \n # A. Data preprocessing:\n1. check datasets for null values: datasets are free of null values.\n\n2. Create labels for the model:\n\nFirst we will extract all the needed labels(GO term ID) from train_terms.tsv file. There are more than 40,000 labels. In order to simplify our model, we will choose the most frequent 500 GO term IDs as labels vor this ver 6 (previously we have selected 1500 GO terms as label).\n\nThen we will fetch the train_terms data for the relevant labels only. and then we will create the The final labels dataframe (label_df),  composed of 500 columns and 142246 entries. \n\n3. Create embeddings for the protein sequences. We will be using pretrained embeddings form protBERT to\n","metadata":{}},{"cell_type":"markdown","source":"## 1. check data set for any missing values","metadata":{}},{"cell_type":"code","source":"#check sequences for null values:\ndf_sequences.isna().sum()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:17.733979Z","iopub.execute_input":"2023-07-12T08:06:17.734386Z","iopub.status.idle":"2023-07-12T08:06:17.837088Z","shell.execute_reply.started":"2023-07-12T08:06:17.734359Z","shell.execute_reply":"2023-07-12T08:06:17.835858Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#check sequences for null values:\ndf_sequences_test.isna().sum()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:18.933319Z","iopub.execute_input":"2023-07-12T08:06:18.933682Z","iopub.status.idle":"2023-07-12T08:06:19.027238Z","shell.execute_reply.started":"2023-07-12T08:06:18.933654Z","shell.execute_reply":"2023-07-12T08:06:19.026046Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#check data for null values:\ntrain_terms.isna().sum()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:19.369378Z","iopub.execute_input":"2023-07-12T08:06:19.369702Z","iopub.status.idle":"2023-07-12T08:06:21.423782Z","shell.execute_reply.started":"2023-07-12T08:06:19.369678Z","shell.execute_reply":"2023-07-12T08:06:21.422621Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_sequences.describe()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:21.42546Z","iopub.execute_input":"2023-07-12T08:06:21.425716Z","iopub.status.idle":"2023-07-12T08:06:21.448986Z","shell.execute_reply.started":"2023-07-12T08:06:21.425693Z","shell.execute_reply":"2023-07-12T08:06:21.44806Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 2. get labels as label_df\nFirst we will extract all the needed labels(GO term ID) from train_terms.tsv file.\n\nThere are more than 40,000 labels. In order to simplify our model, we will choose the most frequent 500 GO term IDs as labels.|","metadata":{}},{"cell_type":"code","source":"# We will now save the first 500 most frequent GO term Ids into a list.\n# Set the limit for label\nnum_of_labels = 500\n\n# Take value counts in descending order and fetch first 1500 `GO term ID` as labels\nlabels = train_terms['term'].value_counts().index[:num_of_labels].tolist()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:40.33887Z","iopub.execute_input":"2023-07-12T08:06:40.339277Z","iopub.status.idle":"2023-07-12T08:06:40.764348Z","shell.execute_reply.started":"2023-07-12T08:06:40.339245Z","shell.execute_reply":"2023-07-12T08:06:40.763294Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Fetch the train_terms data for the relevant labels only\ntrain_terms_updated = train_terms.loc[train_terms['term'].isin(labels)]\ntrain_terms_updated.shape","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:45.319916Z","iopub.execute_input":"2023-07-12T08:06:45.320285Z","iopub.status.idle":"2023-07-12T08:06:45.759247Z","shell.execute_reply.started":"2023-07-12T08:06:45.320258Z","shell.execute_reply":"2023-07-12T08:06:45.758388Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"train_protein_ids = np.load('/kaggle/input/protbert-embeddings-for-cafa5/train_ids.npy')\nprint(train_protein_ids.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:06:49.947877Z","iopub.execute_input":"2023-07-12T08:06:49.948263Z","iopub.status.idle":"2023-07-12T08:06:50.009744Z","shell.execute_reply.started":"2023-07-12T08:06:49.948235Z","shell.execute_reply":"2023-07-12T08:06:50.009086Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Setup progressbar settings.\n# This is strictly for aesthetic.\nbar = progressbar.ProgressBar(maxval=num_of_labels, \\\n    widgets=[progressbar.Bar('=', '[', ']'), ' ', progressbar.Percentage()])\n\n# Create an empty dataframe of required size for storing the labels,\n# i.e, train_size x num_of_labels (142246 x 500)\ntrain_size = train_protein_ids.shape[0] # len(X)\ntrain_labels = np.zeros((train_size ,num_of_labels))\n\n# Convert from numpy to pandas series for better handling\nseries_train_protein_ids = pd.Series(train_protein_ids)\n\n# Loop through each label\nfor i in range(num_of_labels):\n    # For each label, fetch the corresponding train_terms data\n    n_train_terms = train_terms_updated[train_terms_updated['term'] ==  labels[i]]\n    \n    # Fetch all the unique EntryId aka proteins related to the current label(GO term ID)\n    label_related_proteins = n_train_terms['EntryID'].unique()\n    \n    # In the series_train_protein_ids pandas series, if a protein is related\n    # to the current label, then mark it as 1, else 0.\n    # Replace the ith column of train_Y with with that pandas series.\n    train_labels[:,i] =  series_train_protein_ids.isin(label_related_proteins).astype(float)\n    \n    # Progress bar percentage increase\n    bar.update(i+1)\n\n# Notify the end of progress bar \nbar.finish()\n\n# Convert train_Y numpy into pandas dataframe\nlabels_df = pd.DataFrame(data = train_labels, columns = labels)\nprint(labels_df.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:07:15.154451Z","iopub.execute_input":"2023-07-12T08:07:15.154785Z","iopub.status.idle":"2023-07-12T08:09:10.840691Z","shell.execute_reply.started":"2023-07-12T08:07:15.154754Z","shell.execute_reply":"2023-07-12T08:09:10.839826Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The final labels dataframe (label_df) is composed of 500 columns and 142246 entries. We can see all 500 dimensions(results will be truncated since the number of columns is big) of our dataset by printing out the first 5 entries using the following code:","metadata":{}},{"cell_type":"code","source":"labels_df.head()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:26.115901Z","iopub.execute_input":"2023-07-12T08:09:26.116442Z","iopub.status.idle":"2023-07-12T08:09:26.14158Z","shell.execute_reply.started":"2023-07-12T08:09:26.116415Z","shell.execute_reply":"2023-07-12T08:09:26.140699Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Get the sequence embeddings using protBERT","metadata":{}},{"cell_type":"markdown","source":"### What is ProtBert model?\nPretrained model on protein sequences using a masked language modeling (MLM) objective. ProtBert is based on Bert model which pretrained on a large corpus of protein sequences in a self-supervised fashion. This means it was pretrained on the raw protein sequences only, with no humans labelling them in any way (which is why it can use lots of publicly available data) with an automatic process to generate inputs and labels from those protein sequences.\n\nOne important difference between our Bert model and the original Bert version is the way of dealing with sequences as separate documents. This means the Next sentence prediction is not used, as each sequence is treated as a complete document. The masking follows the original Bert training with randomly masks 15% of the amino acids in the input.\n\nAt the end, the feature extracted from this model revealed that the LM-embeddings from unlabeled data (only protein sequences) captured important biophysical properties governing protein shape. This implied learning some of the grammar of the language of life realized in protein sequences.\n\nref: https://huggingface.co/Rostlab/prot_bert\n\nThe precalculated protein embeddings can be found here: https://www.kaggle.com/datasets/henriupton/protbert-embeddings-for-cafa5","metadata":{}},{"cell_type":"code","source":"%%time\ndevice = torch.device('cuda' if torch.cuda.is_available() else 'cpu')\ntokenizer = BertTokenizer.from_pretrained(\"Rostlab/prot_bert\", do_lower_case=False )\nmodel = BertModel.from_pretrained(\"Rostlab/prot_bert\").to(device)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:32.803906Z","iopub.execute_input":"2023-07-12T08:09:32.804467Z","iopub.status.idle":"2023-07-12T08:09:44.008321Z","shell.execute_reply.started":"2023-07-12T08:09:32.804436Z","shell.execute_reply":"2023-07-12T08:09:44.007204Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Lets take a moment to see how amino acids are tokenized by prot-BERT model","metadata":{}},{"cell_type":"code","source":"#get the token embedding for amino acids\ntoken_embedding = {token: model.get_input_embeddings()(torch.tensor(idx)) for token, idx in tokenizer.get_vocab().items()}","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:44.010329Z","iopub.execute_input":"2023-07-12T08:09:44.01061Z","iopub.status.idle":"2023-07-12T08:09:44.02629Z","shell.execute_reply.started":"2023-07-12T08:09:44.010586Z","shell.execute_reply":"2023-07-12T08:09:44.024535Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(len(token_embedding))","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:44.027661Z","iopub.execute_input":"2023-07-12T08:09:44.028042Z","iopub.status.idle":"2023-07-12T08:09:44.032811Z","shell.execute_reply.started":"2023-07-12T08:09:44.027996Z","shell.execute_reply":"2023-07-12T08:09:44.03205Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#print embedding of aminoacid == ALANINE(A)\nprint(token_embedding['A'])","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:44.034626Z","iopub.execute_input":"2023-07-12T08:09:44.035102Z","iopub.status.idle":"2023-07-12T08:09:44.080154Z","shell.execute_reply.started":"2023-07-12T08:09:44.035078Z","shell.execute_reply":"2023-07-12T08:09:44.07946Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(token_embedding)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:44.081087Z","iopub.execute_input":"2023-07-12T08:09:44.081928Z","iopub.status.idle":"2023-07-12T08:09:44.093652Z","shell.execute_reply.started":"2023-07-12T08:09:44.081906Z","shell.execute_reply":"2023-07-12T08:09:44.092765Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We will now load the pre calculated protein embeddings created by HENRI UPTON using the Rost Lab's protBERT protein language model.\n\nThe protein embeddings to be used for training are recorded in train_embeddings.npy and the corresponding protein ids are available in train_ids.npy.\n\nFirst, we will load the protein ids of the protein embeddings in the train dataset contained in train_ids.npy into a numpy array.","metadata":{}},{"cell_type":"markdown","source":"The train_protein_ids array consists of 142246 protein_ids. Let us print out the first 5 entries using the following code:","metadata":{}},{"cell_type":"code","source":"train_protein_ids[:5]","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:44.095388Z","iopub.execute_input":"2023-07-12T08:09:44.096472Z","iopub.status.idle":"2023-07-12T08:09:44.101327Z","shell.execute_reply.started":"2023-07-12T08:09:44.096444Z","shell.execute_reply":"2023-07-12T08:09:44.10067Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"After loading the files as numpy arrays, we will convert them into Pandas dataframe.\n\nEach protein embedding is a vector of length 1024. We create the resulting dataframe such that there are 1024 columns to represent the values in each of the 1024 places in the vector.","metadata":{}},{"cell_type":"code","source":"#get pretrained embeddings:\n\nimport numpy as np\n# Load the .npy file\ndata_embed = np.load('/kaggle/input/protbert-embeddings-for-cafa5/train_embeddings.npy')","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:49.651797Z","iopub.execute_input":"2023-07-12T08:09:49.652146Z","iopub.status.idle":"2023-07-12T08:09:54.624273Z","shell.execute_reply.started":"2023-07-12T08:09:49.65212Z","shell.execute_reply":"2023-07-12T08:09:54.623381Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#get temp df of the embeddings:\ntrain_df = pd.DataFrame(data_embed)\n#del embeddings_list\ntrain_df .columns = ['Feature_' + str(x) for x in train_df.columns]\ntrain_df ","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:55.303879Z","iopub.execute_input":"2023-07-12T08:09:55.304225Z","iopub.status.idle":"2023-07-12T08:09:55.412001Z","shell.execute_reply.started":"2023-07-12T08:09:55.3042Z","shell.execute_reply":"2023-07-12T08:09:55.410106Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"test_embed = np.load('/kaggle/input/protbert-embeddings-for-cafa5/test_embeddings.npy')\n#get temp df for test data embeddings\ntest_df = pd.DataFrame(test_embed)\n#del test_embeddings_list\ntest_df.columns = ['Feature_' + str(x) for x in test_df.columns]\ntest_df","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:09:59.883956Z","iopub.execute_input":"2023-07-12T08:09:59.884321Z","iopub.status.idle":"2023-07-12T08:10:06.990592Z","shell.execute_reply.started":"2023-07-12T08:09:59.884296Z","shell.execute_reply":"2023-07-12T08:10:06.989775Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Model:\nNext, we will use Tensorflow to train a Deep Neural Network with the protein embeddings.","metadata":{}},{"cell_type":"markdown","source":"## 1. define X and y","metadata":{}},{"cell_type":"markdown","source":"## 2. Prepare train-test split","metadata":{}},{"cell_type":"code","source":"# train and validation set\nX = train_df.to_numpy()\ny = labels_df.to_numpy()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:10:08.831743Z","iopub.execute_input":"2023-07-12T08:10:08.832064Z","iopub.status.idle":"2023-07-12T08:10:08.837448Z","shell.execute_reply.started":"2023-07-12T08:10:08.832041Z","shell.execute_reply":"2023-07-12T08:10:08.835803Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(X.shape, y.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:10:10.530326Z","iopub.execute_input":"2023-07-12T08:10:10.530705Z","iopub.status.idle":"2023-07-12T08:10:10.535944Z","shell.execute_reply.started":"2023-07-12T08:10:10.530673Z","shell.execute_reply":"2023-07-12T08:10:10.534973Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"X_train, X_val, y_train, y_val = train_test_split(X, y, test_size=0.2, random_state=42)\n\nprint(X_train.shape, X_val.shape)\nprint(y_train.shape, y_val.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:10:13.470078Z","iopub.execute_input":"2023-07-12T08:10:13.470415Z","iopub.status.idle":"2023-07-12T08:10:13.713812Z","shell.execute_reply.started":"2023-07-12T08:10:13.470389Z","shell.execute_reply":"2023-07-12T08:10:13.712984Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## 3. ANN model","metadata":{}},{"cell_type":"code","source":"INPUT_SHAPE = [X_train.shape[1]]\nBATCH_SIZE = 5120\n\nmodel = tf.keras.Sequential([\n    tf.keras.layers.BatchNormalization(input_shape=INPUT_SHAPE),    \n    tf.keras.layers.Dense(units=512, activation='relu'),\n    tf.keras.layers.Dense(units=512, activation='relu'),\n    tf.keras.layers.Dense(units=512, activation='relu'),\n    tf.keras.layers.Dense(units=num_of_labels,activation='sigmoid')\n])\n\n\n# Compile model\nmodel.compile(\n    optimizer=tf.keras.optimizers.Adam(learning_rate=0.001),\n    loss='binary_crossentropy',\n    metrics=['binary_accuracy', tf.keras.metrics.AUC()],\n)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:11:49.50807Z","iopub.execute_input":"2023-07-12T08:11:49.508444Z","iopub.status.idle":"2023-07-12T08:11:49.768614Z","shell.execute_reply.started":"2023-07-12T08:11:49.508417Z","shell.execute_reply":"2023-07-12T08:11:49.767698Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# callbacks\nearlystopping_cb = tf.keras.callbacks.EarlyStopping(monitor=\"loss\", patience=10)\ncheckpoint_cb = tf.keras.callbacks.ModelCheckpoint(\"/kaggle/working/best_model.pro5\", save_best_only=True)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:11:53.304012Z","iopub.execute_input":"2023-07-12T08:11:53.30442Z","iopub.status.idle":"2023-07-12T08:11:53.310084Z","shell.execute_reply.started":"2023-07-12T08:11:53.304393Z","shell.execute_reply":"2023-07-12T08:11:53.308795Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# train new model\nhistory = model.fit(X_train, y_train, epochs=30,\n                       validation_data=(X_val, y_val),\n                       callbacks=[earlystopping_cb, checkpoint_cb])","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:11:54.827969Z","iopub.execute_input":"2023-07-12T08:11:54.828874Z","iopub.status.idle":"2023-07-12T08:33:45.88842Z","shell.execute_reply.started":"2023-07-12T08:11:54.828838Z","shell.execute_reply":"2023-07-12T08:33:45.887288Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Plot the model's loss and accuracy for each epoch\nA validation accuracy of  0.95\n\n\nThe history.history variable contains the information of each epoch.\n\n","metadata":{}},{"cell_type":"code","source":"history_df = pd.DataFrame(history.history)\nhistory_df.loc[:, ['loss']].plot(title=\"Cross-entropy\")\nhistory_df.loc[:, ['binary_accuracy']].plot(title=\"Accuracy\")","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:34:18.568074Z","iopub.execute_input":"2023-07-12T08:34:18.569301Z","iopub.status.idle":"2023-07-12T08:34:19.027138Z","shell.execute_reply.started":"2023-07-12T08:34:18.569259Z","shell.execute_reply":"2023-07-12T08:34:19.025997Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Submission","metadata":{}},{"cell_type":"markdown","source":"Load features for submission","metadata":{}},{"cell_type":"code","source":"%%time\n# fn = '/kaggle/input/protein-embeddings-1/reduced_embeddings_file.npy'\n# fn = '/kaggle/input/protein-embeddings-1/embed_protbert_train_clip_1200_first_70000_prot.csv'\n# fn = '/kaggle/input/t5embeds/train_embeds.npy'\nfn = '/kaggle/input/protbert-embeddings-for-cafa5/test_embeddings.npy'\nprint(fn)\nX_submit = np.load(fn)\nprint(X_submit.shape)\n# X_submit","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:34:33.931547Z","iopub.execute_input":"2023-07-12T08:34:33.931903Z","iopub.status.idle":"2023-07-12T08:34:34.073969Z","shell.execute_reply.started":"2023-07-12T08:34:33.93188Z","shell.execute_reply":"2023-07-12T08:34:34.072614Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Calculate prediction for submission","metadata":{}},{"cell_type":"code","source":"%%time\nY_submit =  model.predict(X_submit)\nprint(Y_submit.shape)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:18.874487Z","iopub.execute_input":"2023-07-12T08:35:18.874859Z","iopub.status.idle":"2023-07-12T08:35:36.939751Z","shell.execute_reply.started":"2023-07-12T08:35:18.874829Z","shell.execute_reply":"2023-07-12T08:35:36.938873Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time \ndf_finalSubmission = pd.DataFrame(columns = ['Protein Id', 'GO Term Id','Prediction'])","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:36.941357Z","iopub.execute_input":"2023-07-12T08:35:36.941631Z","iopub.status.idle":"2023-07-12T08:35:36.95108Z","shell.execute_reply.started":"2023-07-12T08:35:36.941605Z","shell.execute_reply":"2023-07-12T08:35:36.949733Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Load protein ids for the submission","metadata":{}},{"cell_type":"code","source":"%%time\nfn = '/kaggle/input/protbert-embeddings-for-cafa5/test_ids.npy'\nvec_test_protein_ids = np.load(fn)\nprint(vec_test_protein_ids.shape)\nvec_test_protein_ids","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:36.952509Z","iopub.execute_input":"2023-07-12T08:35:36.952875Z","iopub.status.idle":"2023-07-12T08:35:37.034467Z","shell.execute_reply.started":"2023-07-12T08:35:36.952847Z","shell.execute_reply":"2023-07-12T08:35:37.032954Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\"Melt\" protein ids","metadata":{}},{"cell_type":"code","source":"%%time \nl = []\nfor k in list(vec_test_protein_ids):\n    l += [ k] * Y_submit.shape[1]\nprint(len(l), l[:20])    \n\ndf_finalSubmission['Protein Id'] = l","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:37.037398Z","iopub.execute_input":"2023-07-12T08:35:37.037731Z","iopub.status.idle":"2023-07-12T08:35:52.801442Z","shell.execute_reply.started":"2023-07-12T08:35:37.037702Z","shell.execute_reply":"2023-07-12T08:35:52.799981Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\"Melt\" Labels (Gene ontology terms )","metadata":{}},{"cell_type":"code","source":"df_finalSubmission['GO Term Id'] = labels * Y_submit.shape[0]\n# df_finalSubmission.head(3)","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:52.802586Z","iopub.execute_input":"2023-07-12T08:35:52.802868Z","iopub.status.idle":"2023-07-12T08:35:56.576386Z","shell.execute_reply.started":"2023-07-12T08:35:52.802842Z","shell.execute_reply":"2023-07-12T08:35:56.57564Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Assign predictions","metadata":{}},{"cell_type":"code","source":"df_finalSubmission['Prediction'] = Y_submit.ravel()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:56.577438Z","iopub.execute_input":"2023-07-12T08:35:56.577832Z","iopub.status.idle":"2023-07-12T08:35:56.837387Z","shell.execute_reply.started":"2023-07-12T08:35:56.57781Z","shell.execute_reply":"2023-07-12T08:35:56.836432Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Save","metadata":{}},{"cell_type":"code","source":"%%time\ndf_finalSubmission.to_csv(\"submission.tsv\",header=False, index=False, sep=\"\\t\")","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:35:56.838428Z","iopub.execute_input":"2023-07-12T08:35:56.83867Z","iopub.status.idle":"2023-07-12T08:38:19.825698Z","shell.execute_reply.started":"2023-07-12T08:35:56.838649Z","shell.execute_reply":"2023-07-12T08:38:19.824805Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_finalSubmission.head()","metadata":{"execution":{"iopub.status.busy":"2023-07-12T08:38:19.8268Z","iopub.execute_input":"2023-07-12T08:38:19.827949Z","iopub.status.idle":"2023-07-12T08:38:19.839393Z","shell.execute_reply.started":"2023-07-12T08:38:19.827918Z","shell.execute_reply":"2023-07-12T08:38:19.838379Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}