{"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":"markdown","source":"# Exploratory Data Analysis (EDA) - Gene Ontology (GO) Terms\n","metadata":{}},{"cell_type":"markdown","source":"In this notebook, we will perform an Exploratory Data Analysis (EDA) on the gene ontology terms associated with the sequences in the train dataset. The gene ontology terms provide functional annotations to the sequences, allowing us to gain insights into the biological processes, cellular components, and molecular functions associated with the sequences.\n\n## Dataset Overview\n\nThe analysis will be based on the following files:\n\n- `train_terms.tsv`: This file contains the gene ontology terms associated with each sequence in the train dataset. The terms are represented as a tab-separated values (TSV) file, where each row corresponds to a sequence and the associated gene ontology terms are listed.\n- `go-basic.obo`: This file contains the ontology data in the GO graph structure. The GO graph represents the hierarchical relationships between gene ontology terms, providing information about the parent-child relationships and the overall structure of the gene ontology.\n\n## Goals of the Analysis\n\nThe main goals of this analysis are:\n\n1. Analyze the distribution of GO terms across the sequences to understand the functional annotations present in the dataset.\n2. Identify the most common or prevalent GO terms to gain insights into the dominant functional annotations.\n3. Explore the relationships and dependencies between GO terms using the GO graph structure to uncover the hierarchical organization of functional annotations.\n\n## Approach\n\nTo achieve these goals, we will perform various exploratory data analysis techniques, including data visualization and statistical analysis. We will use Python and relevant libraries such as pandas, matplotlib, and networkx to load, process, and analyze the data.\n\nLet's begin the analysis and delve into the gene ontology terms associated with the sequences in the train dataset!","metadata":{}},{"cell_type":"code","source":"# This Python 3 environment comes with many helpful analytics libraries installed\n# It is defined by the kaggle/python Docker image: https://github.com/kaggle/docker-python\n# For example, here's several helpful packages to load\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport matplotlib.pyplot as plt\n\n# Input data files are available in the read-only \"../input/\" directory\n# For example, running this (by clicking run or pressing Shift+Enter) will list all files under the input directory\n\nimport os\nfor dirname, _, filenames in os.walk('/kaggle/input'):\n    for filename in filenames:\n        print(os.path.join(dirname, filename))\n\n# You can write up to 20GB to the current directory (/kaggle/working/) that gets preserved as output when you create a version using \"Save & Run All\" \n# You can also write temporary files to /kaggle/temp/, but they won't be saved outside of the current session","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2023-05-31T18:23:34.155725Z","iopub.execute_input":"2023-05-31T18:23:34.15617Z","iopub.status.idle":"2023-05-31T18:23:34.180456Z","shell.execute_reply.started":"2023-05-31T18:23:34.156137Z","shell.execute_reply":"2023-05-31T18:23:34.179607Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pandas as pd\n\ntrain_terms_df = pd.read_csv('/kaggle/input/cafa-5-protein-function-prediction/Train/train_terms.tsv', delimiter='\\t')\ngo_terms_df = pd.read_csv('/kaggle/input/cafa-5-protein-function-prediction/Train/train_terms.tsv', delimiter='\\t')\n","metadata":{"execution":{"iopub.status.busy":"2023-05-31T18:31:34.64018Z","iopub.execute_input":"2023-05-31T18:31:34.640611Z","iopub.status.idle":"2023-05-31T18:31:40.563183Z","shell.execute_reply.started":"2023-05-31T18:31:34.640576Z","shell.execute_reply":"2023-05-31T18:31:40.561997Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Perform initial data exploration\nprint(\"Number of rows:\", len(train_terms_df))\nprint(\"Number of columns:\", len(train_terms_df.columns))\nprint(\"Data types:\")\nprint(train_terms_df.dtypes)\nprint(\"Missing values:\")\nprint(train_terms_df.isnull().sum())\n\n# Sample the first few rows of the DataFrame\nprint(\"Sample data:\")\nprint(train_terms_df.head())\n","metadata":{"execution":{"iopub.status.busy":"2023-05-29T18:26:35.98232Z","iopub.execute_input":"2023-05-29T18:26:35.982762Z","iopub.status.idle":"2023-05-29T18:26:38.143063Z","shell.execute_reply.started":"2023-05-29T18:26:35.982729Z","shell.execute_reply":"2023-05-29T18:26:38.141762Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"We can see that the dataset consists of 5,363,863 rows and 3 columns. The columns are labeled as EntryID, term, and aspect, with each column having an object data type.\n\nFurthermore, there are no missing values in any of the columns, as indicated by the count of 0 for each column.\n\nThe sample data shows a subset of the dataset, displaying the EntryID, term, and aspect values for five rows.","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport matplotlib.pyplot as plt\n\n# Read the train_terms.tsv file into a DataFrame\ndata = pd.read_csv('/kaggle/input/cafa-5-protein-function-prediction/Train/train_terms.tsv', sep='\\t')\n\n# Calculate the frequency of different GO terms\ngo_term_counts = data['term'].value_counts()\n\n# Read the go-basic.obo file and extract the mapping of term codes to names\ngo_terms_mapping = {}\nwith open('/kaggle/input/cafa-5-protein-function-prediction/Train/go-basic.obo', 'r') as file:\n    for line in file:\n        if line.startswith('id:'):\n            term_code = line.strip().split(' ')[-1]\n        elif line.startswith('name:'):\n            term_name = line.strip().split(' ')[1:]\n            go_terms_mapping[term_code] = term_name\n\n# Replace the term codes in the index of go_term_counts with the corresponding term names\ngo_term_counts.index = go_term_counts.index.map(go_terms_mapping)\n\n# Plot the distribution of GO terms using a bar plot\nplt.figure(figsize=(10, 6))\ngo_term_counts.head(20).plot(kind='bar')\nplt.xlabel('GO Term')\nplt.ylabel('Frequency')\nplt.title('Distribution of GO Terms (Top 20)')\nplt.xticks(rotation=45, ha='right')\nplt.tight_layout()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-05-31T18:51:20.576746Z","iopub.execute_input":"2023-05-31T18:51:20.577119Z","iopub.status.idle":"2023-05-31T18:51:25.521003Z","shell.execute_reply.started":"2023-05-31T18:51:20.577091Z","shell.execute_reply":"2023-05-31T18:51:25.519917Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"go_term_counts.head(20)","metadata":{"execution":{"iopub.status.busy":"2023-05-29T18:34:07.957337Z","iopub.execute_input":"2023-05-29T18:34:07.95766Z","iopub.status.idle":"2023-05-29T18:34:07.971499Z","shell.execute_reply.started":"2023-05-29T18:34:07.957638Z","shell.execute_reply":"2023-05-29T18:34:07.970009Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"These terms represent different functional categories such as cellular components, biological processes, and molecular functions. Here are some observations and possible next steps based on the prevalent terms:\n\n1. **Cellular Component**: The terms related to cellular components are highly prevalent, indicating that the dataset contains a significant focus on the organization and localization of cellular structures.\n2. **Biological Process**: Similarly, the prevalence of terms related to biological processes suggests that the dataset covers a diverse range of biological activities and functions.\n3. **Molecular Function**: Molecular function terms indicate the types of activities performed by individual molecules, such as binding, catalysis, or receptor activity.\n4. **Intracellular Structures and Organelles**: Terms related to intracellular structures and organelles are prominent, emphasizing the importance of studying the internal organization and functions of cells.\n\nBased on these prevalent terms, some possible next steps could include:\n\n1. **Further Exploration of Specific Functional Categories**: Dive deeper into specific functional categories by analyzing their distributions, exploring relationships between different terms within a category, and investigating their associations with other variables of interest.\n2. **Functional Enrichment Analysis**: Perform functional enrichment analysis to identify overrepresented or enriched functional categories within the dataset. This analysis can provide insights into the functional significance of the dataset and potentially reveal underlying biological processes or pathways.\n3. **Comparative Analysis**: Compare the distribution of GO terms across different subsets or categories within the dataset. This analysis can help identify differences or similarities in functional annotations among different groups or conditions.","metadata":{}},{"cell_type":"code","source":"unique_categories = go_terms_df['term'].unique()\nprint(unique_categories)","metadata":{"execution":{"iopub.status.busy":"2023-05-29T18:43:37.304073Z","iopub.execute_input":"2023-05-29T18:43:37.304443Z","iopub.status.idle":"2023-05-29T18:43:37.695517Z","shell.execute_reply.started":"2023-05-29T18:43:37.304411Z","shell.execute_reply":"2023-05-29T18:43:37.693913Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"category_counts = data['term'].value_counts()\nprint(category_counts)","metadata":{"execution":{"iopub.status.busy":"2023-05-31T18:34:18.385665Z","iopub.execute_input":"2023-05-31T18:34:18.386913Z","iopub.status.idle":"2023-05-31T18:34:19.389888Z","shell.execute_reply.started":"2023-05-31T18:34:18.38686Z","shell.execute_reply":"2023-05-31T18:34:19.38894Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"The unique values in the \"category\" column are represented by the GO term codes, such as 'GO:0008152', 'GO:0034655', etc. Instead of having explicit category names, the dataset contains the GO term codes as categories.\n\nTo investigate the data further and understand the distribution and availability of GO terms, we can follow these steps:\n\n1. Convert GO term codes to category names","metadata":{}},{"cell_type":"markdown","source":"2. Analyze the distribution of GO terms","metadata":{}},{"cell_type":"code","source":"!pip install goatools","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import requests\nfrom goatools import obo_parser\n\n# Download the OBO file if not already downloaded\nobo_file_url = 'http://purl.obolibrary.org/obo/go/go-basic.obo'\nobo_file_path = '/kaggle/input/cafa-5-protein-function-prediction/Train/go-basic.obo'\n\n\n# Load the GO graph from the OBO file\ngo_graph = obo_parser.GODag(obo_file_path)\n\n# Create a dictionary with term names and their IDs\ngo_names_dict = {go_id: go_term.name for go_id, go_term in go_graph.items()}\n","metadata":{"execution":{"iopub.status.busy":"2023-05-31T18:46:28.065233Z","iopub.execute_input":"2023-05-31T18:46:28.065685Z","iopub.status.idle":"2023-05-31T18:46:30.171028Z","shell.execute_reply.started":"2023-05-31T18:46:28.06565Z","shell.execute_reply":"2023-05-31T18:46:30.169945Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from wordcloud import WordCloud\nfrom collections import Counter\n\n# Flatten the list of term names\nterm_names = [term for terms in go_terms_mapping.values() for term in terms]\n\n# Count the frequency of each term name\nterm_counts = Counter(term_names)\n\n# Convert term names and frequencies to a dictionary\nterm_freq_dict = {term: frequency for term, frequency in term_counts.items()}\ndel term_freq_dict['of']\ndel term_freq_dict['to']\ndel term_freq_dict['in']\ndel term_freq_dict['or']\ndel term_freq_dict['by']\ndel term_freq_dict['into']\n\n# Generate word cloud using term names and frequencies\nwordcloud = WordCloud(width=800, height=400, background_color='white').generate_from_frequencies(term_freq_dict)\n\n# Display the word cloud\nplt.figure(figsize=(12, 6))\nplt.imshow(wordcloud, interpolation='bilinear')\nplt.axis('off')\nplt.title('Distribution of GO Terms')\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-05-31T18:56:24.333079Z","iopub.execute_input":"2023-05-31T18:56:24.33359Z","iopub.status.idle":"2023-05-31T18:56:25.424772Z","shell.execute_reply.started":"2023-05-31T18:56:24.333555Z","shell.execute_reply":"2023-05-31T18:56:25.423723Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Get the top 20 words excluding 'of', 'in', and 'to'\ntop_words = [word for word in term_counts.most_common(23) if word[0] not in ['of', 'in', 'to']]\ntop_word_labels = [word[0] for word in top_words[:20]]\ntop_word_counts = [word[1] for word in top_words[:20]]\n\n# Create a bar plot of the top 20 words\nplt.figure(figsize=(12, 6))\nplt.bar(range(len(top_word_labels)), top_word_counts)\nplt.xlabel('Term')\nplt.ylabel('Frequency')\nplt.title('Top 20 Words (Excluding \"of\", \"in\", and \"to\")')\nplt.xticks(range(len(top_word_labels)), top_word_labels, rotation=45, ha='right')\nplt.tight_layout()\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-05-31T19:00:34.957632Z","iopub.execute_input":"2023-05-31T19:00:34.958798Z","iopub.status.idle":"2023-05-31T19:00:35.459251Z","shell.execute_reply.started":"2023-05-31T19:00:34.95875Z","shell.execute_reply":"2023-05-31T19:00:35.45786Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"Lets look at a random Go Term:","metadata":{}},{"cell_type":"code","source":"go_term = 'GO:0005575'  \nmax_child_terms = 10  # Maximum number of child terms to display\n\n# Retrieve the GO term from the GO graph\nterm = go_graph.get(go_term)\n\nif term is not None:\n    # Get the immediate parent terms\n    parent_terms = []\n    for parent_term_id in term.get_all_parents():\n        parent_term = go_graph.get(parent_term_id)\n        if parent_term is not None:\n            parent_terms.append(parent_term.name)\n\n    # Get the immediate child terms\n    child_terms = []\n    for child_term_id in term.get_all_children():\n        child_term = go_graph.get(child_term_id)\n        if child_term is not None:\n            child_terms.append(child_term.name)\n\n    # Print the information\n    print(\"GO Term:\", go_term)\n    print(\"Term Name:\", term.name)\n    print(\"Immediate Parents:\")\n    if parent_terms:\n        print(\"\\n\".join(parent_terms))\n    else:\n        print(\"None\")\n    print(\"Immediate Children:\")\n    if child_terms:\n        if len(child_terms) <= max_child_terms:\n            print(\"\\n\".join(child_terms))\n        else:\n            print(\"\\n\".join(child_terms[:max_child_terms]) + f\"\\n... and {len(child_terms) - max_child_terms} more\")\n    else:\n        print(\"None\")\n    print(\"-\" * 50)\nelse:\n    print(\"GO term not found in the GO graph.\")\n","metadata":{"execution":{"iopub.status.busy":"2023-05-27T20:17:27.371325Z","iopub.execute_input":"2023-05-27T20:17:27.372511Z","iopub.status.idle":"2023-05-27T20:17:27.392899Z","shell.execute_reply.started":"2023-05-27T20:17:27.372465Z","shell.execute_reply":"2023-05-27T20:17:27.391389Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"\nLets categorize the GO terms based on their namespaces (biological process, molecular function, cellular component) and display the top few GO terms within each category. Remember that GO terms represent functional annotations associated with genes or gene products, describing their biological processes, molecular functions, or cellular components.\n\nTo do so, we'll get an overview of the distribution of GO terms across different functional categories. By categorizing and displaying the top GO terms within each category, we can gain insights into the prevalent functional annotations in the dataset and better understand the functional landscape of the proteins.","metadata":{}},{"cell_type":"code","source":"from goatools.obo_parser import GODag\nfrom collections import defaultdict\n\ndef get_go_categories(go_graph):\n    categories = defaultdict(list)\n    for go_id, go_term in go_graph.items():\n        for parent_id in go_term.parents:\n            if parent_id in go_graph:\n                parent_term = go_graph[parent_id]\n                if parent_term.namespace != go_term.namespace:\n                    categories[parent_id].append(go_id)\n    return categories\n\n\n# Load the GO graph from the OBO file\n#go_graph = GODag('go-basic.obo')\n\n# Derive categories from the GO graph\ncategories = get_go_categories(go_graph)\n","metadata":{"execution":{"iopub.status.busy":"2023-05-27T20:22:17.676391Z","iopub.execute_input":"2023-05-27T20:22:17.676767Z","iopub.status.idle":"2023-05-27T20:22:17.72131Z","shell.execute_reply.started":"2023-05-27T20:22:17.676737Z","shell.execute_reply":"2023-05-27T20:22:17.71929Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from collections import defaultdict\n\n# Derive categories from the GO graph\ncategories = defaultdict(list)\n\nfor go_id, go_term in go_graph.items():\n    if go_term.namespace == 'biological_process':\n        categories['Biological Process'].append(go_term.name)\n    elif go_term.namespace == 'molecular_function':\n        categories['Molecular Function'].append(go_term.name)\n    elif go_term.namespace == 'cellular_component':\n        categories['Cellular Component'].append(go_term.name)\n\n# Define the number of top GO terms to display per category\ntop_n = 10\n\n# Sort the terms within each category by prevalence in the go_term dataframe\nfor category, go_terms in categories.items():\n    go_terms_sorted = sorted(go_terms, key=lambda x: go_term_counts.get(x, 0), reverse=True)\n    print(f\"Category: {category}\")\n    print(\"Top GO Terms:\")\n    for go_term in go_terms_sorted[:top_n]:\n        print(go_term)\n    print(\"-\" * 50)\n\n","metadata":{"execution":{"iopub.status.busy":"2023-05-27T20:30:54.333263Z","iopub.execute_input":"2023-05-27T20:30:54.333665Z","iopub.status.idle":"2023-05-27T20:30:54.568279Z","shell.execute_reply.started":"2023-05-27T20:30:54.333635Z","shell.execute_reply":"2023-05-27T20:30:54.566442Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}