{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"pygments_lexer":"ipython3","nbconvert_exporter":"python","version":"3.6.4","file_extension":".py","codemirror_mode":{"name":"ipython","version":3},"name":"python","mimetype":"text/x-python"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":41875,"databundleVersionId":5521661,"sourceType":"competition"},{"sourceId":5599516,"sourceType":"datasetVersion","datasetId":3221195},{"sourceId":6046570,"sourceType":"datasetVersion","datasetId":3458902},{"sourceId":6046672,"sourceType":"datasetVersion","datasetId":3391266},{"sourceId":6061116,"sourceType":"datasetVersion","datasetId":3468459},{"sourceId":6063292,"sourceType":"datasetVersion","datasetId":3469947},{"sourceId":6064056,"sourceType":"datasetVersion","datasetId":3470430},{"sourceId":6125456,"sourceType":"datasetVersion","datasetId":3511563},{"sourceId":6129570,"sourceType":"datasetVersion","datasetId":3514271},{"sourceId":6311022,"sourceType":"datasetVersion","datasetId":3630996},{"sourceId":6319967,"sourceType":"datasetVersion","datasetId":3636951},{"sourceId":6326388,"sourceType":"datasetVersion","datasetId":3641107},{"sourceId":6331075,"sourceType":"datasetVersion","datasetId":3644155},{"sourceId":6336053,"sourceType":"datasetVersion","datasetId":3647520},{"sourceId":6336913,"sourceType":"datasetVersion","datasetId":3648012},{"sourceId":6342556,"sourceType":"datasetVersion","datasetId":3651794},{"sourceId":7365880,"sourceType":"datasetVersion","datasetId":4279154},{"sourceId":7379084,"sourceType":"datasetVersion","datasetId":4288208},{"sourceId":7384998,"sourceType":"datasetVersion","datasetId":4292368},{"sourceId":7386015,"sourceType":"datasetVersion","datasetId":4293072},{"sourceId":133778155,"sourceType":"kernelVersion"}],"dockerImageVersionId":30458,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"raw","source":"","metadata":{}},{"cell_type":"markdown","source":"# Combining datasets\n\n### Here we work towards a method to combine datasets and outputs predictions\n\n#### We will use: SPROF Predictions, ProFun (https://github.com/SamusRam/ProFun) Predictions and QuickGo Annotations","metadata":{}},{"cell_type":"code","source":"#Create a dictionary to assign each go term to the roots (CCO, MFO, BPO)\n\nimport re\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport joblib\nimport pickle\nfrom tqdm import tqdm\nfrom Bio import SeqIO\nimport gc\nimport matplotlib.pyplot as plt\n\ndef extract_go_terms_and_branches(file_path):\n    with open(file_path, 'r') as file:\n        content = file.read()\n        # Match each stanza with [Term] in the OBO file\n        stanzas = re.findall(r'\\[Term\\][\\s\\S]*?(?=\\n\\[|$)', content)\n\n    go_terms_dict = {}\n    for stanza in stanzas:\n        # Extract the GO term ID\n        go_id = re.search(r'^id: (GO:\\d+)', stanza, re.MULTILINE)\n        if go_id:\n            go_id = go_id.group(1)\n\n        # Extract the namespace (branch)\n        namespace = re.search(r'^namespace: (\\w+)', stanza, re.MULTILINE)\n        if namespace:\n            namespace = namespace.group(1)\n\n        if go_id and namespace:\n            # Map the branch abbreviation to the corresponding BPO, CCO, or MFO\n            branch_abbr = {'biological_process': 'BPO', 'cellular_component': 'CCO', 'molecular_function': 'MFO'}\n            go_terms_dict[go_id] = branch_abbr[namespace]\n\n    return go_terms_dict\n\nfile_path = '/kaggle/input/cafa-5-protein-function-prediction/Train/go-basic.obo'\ngo_terms_dict = extract_go_terms_and_branches(file_path)","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:50:42.132214Z","iopub.execute_input":"2023-07-01T19:50:42.132925Z","iopub.status.idle":"2023-07-01T19:50:44.123627Z","shell.execute_reply.started":"2023-07-01T19:50:42.132887Z","shell.execute_reply":"2023-07-01T19:50:44.122414Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# ","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:50:44.126264Z","iopub.execute_input":"2023-07-01T19:50:44.126629Z","iopub.status.idle":"2023-07-01T19:50:44.131388Z","shell.execute_reply.started":"2023-07-01T19:50:44.12659Z","shell.execute_reply":"2023-07-01T19:50:44.129627Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Define a class to manage predictions for proteins.\n# The class keeps track of the highest score for each GO (Gene Ontology) term prediction.\n# Note: This assumes scores are comparable, which might not be the case.\n# A ranking-based selection could be more suitable.\n# Each branch outputs a maximum of 35 predictions for each protein after sorting predictions from highest to lowest.\n# There is an option to add a bonus to the score if the term is predicted by multiple methods.\n\nclass ProteinPredictions:\n    # Initialize an empty dictionary to store the predictions\n    def __init__(self):\n        self.predictions = {}\n\n    # Add a prediction to the storage, with optional bonus\n    # Arguments:\n    #   - protein: Identifier for the protein\n    #   - go_term: GO term that is being predicted\n    #   - score: Confidence score of the prediction\n    #   - branch: Branch of the Gene Ontology (e.g., 'CCO', 'MFO', 'BPO')\n    #   - bonus: Optional bonus to be added to the score\n    def add_prediction(self, protein, go_term, score, branch, bonus=1, adjustment=1, power=3):\n        # If the protein is not already in the storage, initialize its structure\n        if protein not in self.predictions:\n            self.predictions[protein] = {'CCO': {}, 'MFO': {}, 'BPO': {}}\n\n        # Convert the score to a float for comparison and calculation\n        score = float(score)\n\n        # If this GO term has already been predicted for this protein and branch,\n        # add the bonus to the score. Keep the highest score.\n        if go_term in self.predictions[protein][branch]:\n            self.predictions[protein][branch][go_term] *= 1+(score**power)*bonus\n            self.predictions[protein][branch][go_term] += score*adjustment\n\n        # If this GO term has not been predicted yet, store it with the score\n        else:\n            self.predictions[protein][branch][go_term] = max(score*adjustment, 0)\n\n        # Ensure that the score does not exceed 1\n        if self.predictions[protein][branch][go_term] > 1:\n            self.predictions[protein][branch][go_term] = 1\n            \n    # Get a list of all scores in the predictions\n    def get_scores(self):\n        scores = []\n        for protein, branches in self.predictions.items():\n            for branch, go_terms in branches.items():\n                scores.extend(go_terms.values())\n        return scores\n         \n    def plot_predictions(self):\n        scores = self.get_scores()\n        plt.hist(scores, bins=30, edgecolor='black')\n        plt.title('Distribution of Prediction Scores')\n        plt.xlabel('Score')\n        plt.ylabel('Frequency')\n        plt.show()\n        \n    # Export the stored predictions to a file\n    # Arguments:\n    #   - output_file: File name for the exported predictions\n    #   - top: Number of top predictions to export for each protein and branch\n\n    def get_predictions(self, output_file='submission.tsv', top=60):\n        with open(output_file, 'w') as f:\n            for protein, branches in self.predictions.items():\n                # Only consider BPO branch for output\n                bpo_terms = branches.get('MFO', {})\n                top_bpo_terms = sorted(bpo_terms.items(), key=lambda x: x[1], reverse=True)[:top]\n                for go_term, score in top_bpo_terms:\n                    f.write(f\"{protein}\\t{go_term}\\t{score:.3f}\\n\")\n","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:50:44.133089Z","iopub.execute_input":"2023-07-01T19:50:44.134252Z","iopub.status.idle":"2023-07-01T19:50:44.153217Z","shell.execute_reply.started":"2023-07-01T19:50:44.134211Z","shell.execute_reply":"2023-07-01T19:50:44.152265Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"protein_predictions = ProteinPredictions()","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:50:44.15576Z","iopub.execute_input":"2023-07-01T19:50:44.156313Z","iopub.status.idle":"2023-07-01T19:50:44.174366Z","shell.execute_reply.started":"2023-07-01T19:50:44.156281Z","shell.execute_reply":"2023-07-01T19:50:44.173292Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for l in tqdm(open('/kaggle/input/f02-sub/f02_BPO_m25.tsv')):\n    item_list = l.split('\\t')\n    temp_id = item_list[0]\n    go=item_list[1]\n    score = float(item_list[2].strip())\n    if go in go_terms_dict:\n        root = go_terms_dict[go]\n        #branch = item_list[3].strip()\n        protein_predictions.add_prediction(temp_id, go, score, root)\nprotein_predictions.plot_predictions()\n","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:50:44.17605Z","iopub.execute_input":"2023-07-01T19:50:44.176875Z","iopub.status.idle":"2023-07-01T19:52:55.822837Z","shell.execute_reply.started":"2023-07-01T19:50:44.176837Z","shell.execute_reply":"2023-07-01T19:52:55.821855Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for l in tqdm(open('/kaggle/input/057-merged-data/submission.tsv')):\n    item_list = l.split('\\t')\n    if len(item_list) == 3:\n        temp_id = item_list[0]\n        go=item_list[1]\n        score = float(item_list[2].strip())\n        if go in go_terms_dict:\n            root = go_terms_dict[go]\n            #branch = item_list[3].strip()\n            protein_predictions.add_prediction(temp_id, go, score, root, bonus=1, adjustment=0.8, power=3)\nprotein_predictions.plot_predictions()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for l in tqdm(open('/kaggle/input/sprof-quickgo/quickgo_sub.tsv')):\n    item_list = l.split('\\t')\n    temp_id = item_list[1]\n    go=item_list[2].strip()\n    score = float(1)\n    if go in go_terms_dict:\n        root = go_terms_dict[go]\n        #branch = item_list[3].strip()\n        protein_predictions.add_prediction(temp_id, go, score, root, bonus=0.1, adjustment=0.1, power=2)\nprotein_predictions.plot_predictions()","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:52:55.824374Z","iopub.execute_input":"2023-07-01T19:52:55.824729Z","iopub.status.idle":"2023-07-01T19:54:08.392135Z","shell.execute_reply.started":"2023-07-01T19:52:55.824694Z","shell.execute_reply":"2023-07-01T19:54:08.390854Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# for l in tqdm(open('/kaggle/input/sprof-quickgo/sprof_submission.tsv')):\n#     item_list = l.split('\\t')\n#     temp_id = item_list[0]\n#     go=item_list[1]\n#     score = float(item_list[2].strip())\n#     if go in go_terms_dict:\n#         root = go_terms_dict[go]\n#         #branch = item_list[3].strip()\n#         protein_predictions.add_prediction(temp_id, go, score, root, bonus=0.1, adjustment=0.01, power=3)\n# protein_predictions.plot_predictions()","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:54:08.393372Z","iopub.execute_input":"2023-07-01T19:54:08.393659Z","iopub.status.idle":"2023-07-01T19:55:23.704157Z","shell.execute_reply.started":"2023-07-01T19:54:08.393629Z","shell.execute_reply":"2023-07-01T19:55:23.703371Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# for l in tqdm(open('/kaggle/input/deepgoplusoutput-v1/submission.tsv')):\n#     item_list = l.split('\\t')\n#     if len(item_list) == 3:\n#         temp_id = item_list[0]\n#         go=item_list[1]\n#         score = float(item_list[2].strip())\n#         if go in go_terms_dict:\n#             root = go_terms_dict[go]\n#             #branch = item_list[3].strip()\n#             protein_predictions.add_prediction(temp_id, go, score, root, bonus=0.1, adjustment=0.05, power=3)\n# protein_predictions.plot_predictions()","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"item_list","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:58:43.350812Z","iopub.execute_input":"2023-07-01T19:58:43.351218Z","iopub.status.idle":"2023-07-01T19:58:43.358574Z","shell.execute_reply.started":"2023-07-01T19:58:43.351186Z","shell.execute_reply":"2023-07-01T19:58:43.357228Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"protein_predictions.get_predictions()","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:55:23.784109Z","iopub.status.idle":"2023-07-01T19:55:23.784412Z","shell.execute_reply.started":"2023-07-01T19:55:23.784261Z","shell.execute_reply":"2023-07-01T19:55:23.784277Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!head -n 200 'submission.tsv'","metadata":{"execution":{"iopub.status.busy":"2023-07-01T19:55:23.785487Z","iopub.status.idle":"2023-07-01T19:55:23.785784Z","shell.execute_reply.started":"2023-07-01T19:55:23.785635Z","shell.execute_reply":"2023-07-01T19:55:23.785651Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}