{"metadata":{"kernelspec":{"name":"python3","display_name":"Python 3","language":"python"},"language_info":{"name":"python","version":"3.11.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"nvidiaTeslaT4","dataSources":[{"sourceId":116062,"databundleVersionId":14875579,"sourceType":"competition"}],"dockerImageVersionId":31193,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"## 1. Import Libraries","metadata":{}},{"cell_type":"code","source":"import pandas as pd\nimport numpy as np\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nimport os\nfrom collections import Counter\nfrom pathlib import Path\n\n# Set plotting style\nsns.set_style('whitegrid')\nplt.rcParams['figure.figsize'] = (14, 8)\nplt.rcParams['font.size'] = 11\n\nprint(\"✓ Libraries imported successfully\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:14.415492Z","iopub.execute_input":"2025-12-17T10:29:14.415677Z","iopub.status.idle":"2025-12-17T10:29:19.356879Z","shell.execute_reply.started":"2025-12-17T10:29:14.415661Z","shell.execute_reply":"2025-12-17T10:29:19.356068Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 2. Load Data from CAFA06","metadata":{}},{"cell_type":"code","source":"# Set data paths\nDATA_DIR = '/kaggle/input/cafa-6-protein-function-prediction'\nTRAIN_DIR = os.path.join(DATA_DIR, 'Train')\nTEST_DIR = os.path.join(DATA_DIR, 'Test')\n\n# Load data\nprint(\"Loading data files...\")\ntrain_terms = pd.read_csv(os.path.join(TRAIN_DIR, 'train_terms.tsv'), sep='\\t')\ntrain_taxonomy = pd.read_csv(os.path.join(TRAIN_DIR, 'train_taxonomy.tsv'), sep='\\t', header=None, names=['EntryID', 'TaxonomyID'])\nia_scores = pd.read_csv(os.path.join(DATA_DIR, 'IA.tsv'), sep='\\t', header=None, names=['GO_term', 'IA_score'])\ntest_taxonomy = pd.read_csv(os.path.join(TEST_DIR, 'testsuperset-taxon-list.tsv'), sep='\\t')\n\n# Parse FASTA files\ndef parse_fasta(filepath):\n    \"\"\"Parse FASTA file\"\"\"\n    sequences = {}\n    current_id = None\n    current_seq = []\n    \n    with open(filepath, 'r') as f:\n        for line in f:\n            line = line.strip()\n            if line.startswith('>'):\n                if current_id is not None:\n                    sequences[current_id] = ''.join(current_seq)\n                current_id = line[1:].split()[0]\n                current_seq = []\n            else:\n                current_seq.append(line)\n        if current_id is not None:\n            sequences[current_id] = ''.join(current_seq)\n    \n    return sequences\n\nprint(\"Parsing FASTA files...\")\ntrain_sequences = parse_fasta(os.path.join(TRAIN_DIR, 'train_sequences.fasta'))\ntest_sequences = parse_fasta(os.path.join(TEST_DIR, 'testsuperset.fasta'))\nprint(f\"  - Train proteins: {len(train_sequences):,}\")\nprint(f\"  - Test proteins: {len(test_sequences):,}\")\nprint(f\"  - Train annotations: {len(train_terms):,}\")\nprint(f\"  - Unique GO terms: {train_terms['term'].nunique():,}\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:19.358258Z","iopub.execute_input":"2025-12-17T10:29:19.358534Z","iopub.status.idle":"2025-12-17T10:29:22.506827Z","shell.execute_reply.started":"2025-12-17T10:29:19.358516Z","shell.execute_reply":"2025-12-17T10:29:22.506051Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 3. Calculate Basic Statistics","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SECTION 1: SỐ LƯỢNG VÀ CẤU TRÚC DỮ LIỆU\")\nprint(\"=\"*80)\n\n# Basic counts\nn_train_proteins = len(train_sequences)\nn_test_proteins = len(test_sequences)\nn_annotations = len(train_terms)\nn_unique_go_terms = train_terms['term'].nunique()\nn_unique_species = train_taxonomy['TaxonomyID'].nunique()\nn_unique_proteins = train_terms['EntryID'].nunique()\n\n# File sizes\ntrain_fasta_size = os.path.getsize(os.path.join(TRAIN_DIR, 'train_sequences.fasta')) / (1024**2)  # MB\ntrain_terms_size = os.path.getsize(os.path.join(TRAIN_DIR, 'train_terms.tsv')) / (1024**2)  # MB\n\n# Create results dataframe\nresults_df = pd.DataFrame({\n    'Metric': [\n        'Total training proteins',\n        'Total test proteins',\n        'Unique proteins (in train)',\n        'Total annotations (train_terms)',\n        'Unique GO terms',\n        'Unique species (TaxonomyID)',\n        'File size (train_sequences.fasta)',\n        'File size (train_terms.tsv)'\n    ],\n    'Value': [\n        f\"{n_train_proteins:,}\",\n        f\"{n_test_proteins:,}\",\n        f\"{n_unique_proteins:,}\",\n        f\"{n_annotations:,}\",\n        f\"{n_unique_go_terms:,}\",\n        f\"{n_unique_species}\",\n        f\"{train_fasta_size:.1f} MB\",\n        f\"{train_terms_size:.1f} MB\"\n    ]\n})\n\nprint(\"\\nBasic Data Statistics:\")\nprint(results_df.to_string(index=False))\nprint()\n\n# Store for later\nstats_dict = {\n    'n_train_proteins': n_train_proteins,\n    'n_test_proteins': n_test_proteins,\n    'n_annotations': n_annotations,\n    'n_unique_go_terms': n_unique_go_terms\n}","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:22.511626Z","iopub.execute_input":"2025-12-17T10:29:22.511846Z","iopub.status.idle":"2025-12-17T10:29:22.60613Z","shell.execute_reply.started":"2025-12-17T10:29:22.511828Z","shell.execute_reply":"2025-12-17T10:29:22.605551Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 4. Aspect Distribution (C, F, P)","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SECTION 2: PHÂN BỔ THEO ASPECT (GO ANNOTATION TYPE)\")\nprint(\"=\"*80)\n\n# Count by aspect\naspect_counts = train_terms['aspect'].value_counts().sort_index()\naspect_labels = {'C': 'Cellular Component', 'F': 'Molecular Function', 'P': 'Biological Process'}\n\naspect_df = pd.DataFrame({\n    'Aspect': aspect_counts.index,\n    'Description': [aspect_labels[a] for a in aspect_counts.index],\n    'Count': aspect_counts.values,\n    '% of Total': (aspect_counts.values / aspect_counts.sum() * 100).round(1)\n})\n\nprint(\"\\nAspect Distribution:\")\nprint(aspect_df.to_string(index=False))\n\n# Add total row\nprint(f\"\\nTOTAL: {aspect_counts.sum():,} (100.0%)\")\n\n# Visualize\nfig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5))\n\n# Bar chart\ncolors = ['#FF6B6B', '#4ECDC4', '#45B7D1']\nax1.bar(aspect_counts.index, aspect_counts.values, color=colors, alpha=0.7, edgecolor='black', linewidth=2)\nax1.set_xlabel('Aspect', fontsize=12, fontweight='bold')\nax1.set_ylabel('Count', fontsize=12, fontweight='bold')\nax1.set_title('GO Annotations by Aspect', fontsize=13, fontweight='bold')\nax1.set_xticklabels([aspect_labels[a] for a in aspect_counts.index], rotation=15, ha='right')\n\n# Add value labels on bars\nfor i, (idx, val) in enumerate(zip(aspect_counts.index, aspect_counts.values)):\n    pct = val / aspect_counts.sum() * 100\n    ax1.text(i, val, f'{val:,}\\n({pct:.1f}%)', ha='center', va='bottom', fontweight='bold')\n\n# Pie chart\nax2.pie(aspect_counts.values, labels=[aspect_labels[a] for a in aspect_counts.index], \n        autopct='%1.1f%%', colors=colors, startangle=90, textprops={'fontsize': 11, 'fontweight': 'bold'})\nax2.set_title('Aspect Distribution (%)', fontsize=13, fontweight='bold')\n\nplt.tight_layout()\nplt.show()\n\nprint(f\"\\n✓ Aspect distribution visualized\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:22.606813Z","iopub.execute_input":"2025-12-17T10:29:22.607311Z","iopub.status.idle":"2025-12-17T10:29:23.095924Z","shell.execute_reply.started":"2025-12-17T10:29:22.607292Z","shell.execute_reply":"2025-12-17T10:29:23.095215Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 5. GO Term Distribution","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SECTION 3: ĐỘ PHÂN BỐ GO TERMS\")\nprint(\"=\"*80)\n\n# Count terms\ngo_counts = train_terms['term'].value_counts()\n\n# Statistics\ngo_stats_df = pd.DataFrame({\n    'Metric': [\n        'GO terms (unique)',\n        'Most common term',\n        'Most common term count',\n        'Least common term',\n        'Least common term count',\n        'Median annotations/GO term',\n        'Mean annotations/GO term',\n        'Std Dev annotations/GO term'\n    ],\n    'Value': [\n        f\"{go_counts.shape[0]:,}\",\n        f\"{go_counts.index[0]}\",\n        f\"{go_counts.iloc[0]:,} ({go_counts.iloc[0]/len(train_terms)*100:.1f}%)\",\n        f\"{go_counts.index[-1]}\",\n        f\"{go_counts.iloc[-1]}\",\n        f\"{go_counts.median():.0f}\",\n        f\"{go_counts.mean():.1f}\",\n        f\"{go_counts.std():.1f}\"\n    ]\n})\n\nprint(\"\\nGO Term Distribution Statistics:\")\nprint(go_stats_df.to_string(index=False))\n\n# Visualize distribution\nfig, axes = plt.subplots(2, 2, figsize=(14, 10))\n\n# Top 20 most common GO terms\nax = axes[0, 0]\ntop_20 = go_counts.head(20)\nax.barh(range(len(top_20)), top_20.values, color='steelblue', alpha=0.7, edgecolor='black')\nax.set_yticks(range(len(top_20)))\nax.set_yticklabels(top_20.index, fontsize=9)\nax.set_xlabel('Count', fontweight='bold')\nax.set_title('Top 20 Most Common GO Terms', fontweight='bold')\nax.invert_yaxis()\n\n# Distribution histogram\nax = axes[0, 1]\nax.hist(go_counts.values, bins=50, color='coral', alpha=0.7, edgecolor='black')\nax.axvline(go_counts.mean(), color='red', linestyle='--', linewidth=2, label=f'Mean: {go_counts.mean():.1f}')\nax.axvline(go_counts.median(), color='green', linestyle='--', linewidth=2, label=f'Median: {go_counts.median():.0f}')\nax.set_xlabel('Annotations per GO Term', fontweight='bold')\nax.set_ylabel('Frequency', fontweight='bold')\nax.set_title('Distribution of GO Term Frequencies', fontweight='bold')\nax.legend()\nax.set_yscale('log')\n\n# CDF plot\nax = axes[1, 0]\nsorted_counts = sorted(go_counts.values, reverse=True)\ncumulative = np.cumsum(sorted_counts) / np.sum(sorted_counts) * 100\nax.plot(range(len(cumulative)), cumulative, linewidth=2, color='darkblue')\nax.axhline(80, color='red', linestyle='--', alpha=0.5, label='80% threshold')\nax.axhline(90, color='orange', linestyle='--', alpha=0.5, label='90% threshold')\nax.set_xlabel('GO Terms (ranked)', fontweight='bold')\nax.set_ylabel('Cumulative % of Annotations', fontweight='bold')\nax.set_title('Cumulative Distribution of Annotations', fontweight='bold')\nax.legend()\nax.grid(True, alpha=0.3)\n\n# Box plot\nax = axes[1, 1]\nbp = ax.boxplot([go_counts.values], labels=['GO Terms'], patch_artist=True)\nbp['boxes'][0].set_facecolor('lightblue')\nax.set_ylabel('Annotations per GO Term', fontweight='bold')\nax.set_title('Box Plot of GO Term Frequencies', fontweight='bold')\nax.grid(True, alpha=0.3, axis='y')\n\nplt.tight_layout()\nplt.show()\n\nprint(f\"\\n✓ GO term distribution visualized\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:23.096787Z","iopub.execute_input":"2025-12-17T10:29:23.097047Z","iopub.status.idle":"2025-12-17T10:29:24.512031Z","shell.execute_reply.started":"2025-12-17T10:29:23.097029Z","shell.execute_reply":"2025-12-17T10:29:24.511251Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 6. Protein Annotation Distribution","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SECTION 4: ĐỘ PHÂN BỐ THEO PROTEIN\")\nprint(\"=\"*80)\n\n# Count annotations per protein\nprotein_counts = train_terms['EntryID'].value_counts()\n\n# Binning\nbins = [0, 1, 5, 10, float('inf')]\nlabels = ['1 term', '2-5 terms', '6-10 terms', '>10 terms']\nprotein_bins = pd.cut(protein_counts.values, bins=bins, labels=labels, right=True)\nprotein_bin_counts = protein_bins.value_counts().sort_index()\n\n# Statistics\nprotein_stats_df = pd.DataFrame({\n    'Metric': [\n        'Protein with 1 GO term',\n        'Protein with 2–5 terms',\n        'Protein with 6–10 terms',\n        'Protein with >10 terms',\n        '',\n        'Max annotations/protein',\n        'Mean annotations/protein',\n        'Median annotations/protein',\n        'Std Dev annotations/protein'\n    ],\n    'Value': [\n        f\"{protein_bin_counts[0]:,} ({protein_bin_counts[0]/len(protein_counts)*100:.1f}%)\",\n        f\"{protein_bin_counts[1]:,} ({protein_bin_counts[1]/len(protein_counts)*100:.1f}%)\",\n        f\"{protein_bin_counts[2]:,} ({protein_bin_counts[2]/len(protein_counts)*100:.1f}%)\",\n        f\"{protein_bin_counts[3]:,} ({protein_bin_counts[3]/len(protein_counts)*100:.1f}%)\",\n        '',\n        f\"{protein_counts.max()}\",\n        f\"{protein_counts.mean():.1f}\",\n        f\"{protein_counts.median():.0f}\",\n        f\"{protein_counts.std():.1f}\"\n    ]\n})\n\nprint(\"\\nProtein Annotation Distribution:\")\nprint(protein_stats_df.to_string(index=False))\n\n# Visualize\nfig, axes = plt.subplots(2, 2, figsize=(14, 10))\n\n# Bar chart for binned distribution\nax = axes[0, 0]\ncolors_bin = ['#FF6B6B', '#4ECDC4', '#45B7D1', '#96CEB4']\nbars = ax.bar(range(len(protein_bin_counts)), protein_bin_counts.values, color=colors_bin, alpha=0.7, edgecolor='black', linewidth=2)\nax.set_xticks(range(len(protein_bin_counts)))\nax.set_xticklabels(protein_bin_counts.index, rotation=15, ha='right')\nax.set_ylabel('Number of Proteins', fontweight='bold')\nax.set_title('Distribution of Annotations per Protein', fontweight='bold')\n\n# Add labels\nfor i, (bar, val) in enumerate(zip(bars, protein_bin_counts.values)):\n    pct = val / protein_bin_counts.sum() * 100\n    ax.text(bar.get_x() + bar.get_width()/2, bar.get_height(), \n           f'{val:,}\\n({pct:.1f}%)', ha='center', va='bottom', fontweight='bold', fontsize=9)\n\n# Histogram\nax = axes[0, 1]\nax.hist(protein_counts.values, bins=50, color='skyblue', alpha=0.7, edgecolor='black')\nax.axvline(protein_counts.mean(), color='red', linestyle='--', linewidth=2, label=f'Mean: {protein_counts.mean():.1f}')\nax.axvline(protein_counts.median(), color='green', linestyle='--', linewidth=2, label=f'Median: {protein_counts.median():.0f}')\nax.set_xlabel('Annotations per Protein', fontweight='bold')\nax.set_ylabel('Frequency', fontweight='bold')\nax.set_title('Histogram of Protein Annotations', fontweight='bold')\nax.set_yscale('log')\nax.legend()\n\n# Top 20 proteins\nax = axes[1, 0]\ntop_proteins = protein_counts.head(20)\nax.barh(range(len(top_proteins)), top_proteins.values, color='mediumpurple', alpha=0.7, edgecolor='black')\nax.set_yticks(range(len(top_proteins)))\nax.set_yticklabels(top_proteins.index, fontsize=8)\nax.set_xlabel('Number of GO Terms', fontweight='bold')\nax.set_title('Top 20 Proteins with Most GO Annotations', fontweight='bold')\nax.invert_yaxis()\n\n# Box plot\nax = axes[1, 1]\nbp = ax.boxplot([protein_counts.values], labels=['Proteins'], patch_artist=True)\nbp['boxes'][0].set_facecolor('lightgreen')\nax.set_ylabel('GO Terms per Protein', fontweight='bold')\nax.set_title('Box Plot of Protein Annotations', fontweight='bold')\nax.grid(True, alpha=0.3, axis='y')\n\nplt.tight_layout()\nplt.show()\n\nprint(f\"\\n✓ Protein distribution visualized\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:24.513607Z","iopub.execute_input":"2025-12-17T10:29:24.514096Z","iopub.status.idle":"2025-12-17T10:29:25.850252Z","shell.execute_reply.started":"2025-12-17T10:29:24.514077Z","shell.execute_reply":"2025-12-17T10:29:25.849506Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 7. IA Scores Distribution (Rarity of GO Terms)","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SECTION 5: IA SCORES – RARITY OF GO TERMS\")\nprint(\"=\"*80)\n\n# Define bins for IA scores\nia_bins = [0, 1, 2, 3, 4, 6]\nia_labels = ['0.0–1.0', '1.0–2.0', '2.0–3.0', '3.0–4.0', '4.0–6.0']\nia_type_labels = ['Very Common', 'Common', 'Semi-Rare', 'Rare', 'Very Rare']\n\nia_scores['IA_category'] = pd.cut(ia_scores['IA_score'], bins=ia_bins, labels=ia_labels, right=False)\nia_category_counts = ia_scores['IA_category'].value_counts().sort_index()\n\n# Statistics\nia_stats_df = pd.DataFrame({\n    'IA Score Range': ia_labels,\n    'Count': [ia_category_counts.get(label, 0) for label in ia_labels],\n    '%': [(ia_category_counts.get(label, 0) / len(ia_scores) * 100) for label in ia_labels],\n    'Type': ia_type_labels\n})\n\nprint(\"\\nIA Score Distribution:\")\nia_stats_display = ia_stats_df.copy()\nia_stats_display['%'] = ia_stats_display['%'].round(1)\nprint(ia_stats_display.to_string(index=False))\n\nprint(f\"\\nOverall IA Score Statistics:\")\nprint(f\"  - Min IA Score: {ia_scores['IA_score'].min():.4f}\")\nprint(f\"  - Max IA Score: {ia_scores['IA_score'].max():.4f}\")\nprint(f\"  - Mean IA Score: {ia_scores['IA_score'].mean():.4f}\")\nprint(f\"  - Median IA Score: {ia_scores['IA_score'].median():.4f}\")\nprint(f\"  - Std Dev: {ia_scores['IA_score'].std():.4f}\")\n\n# Visualize\nfig, axes = plt.subplots(2, 2, figsize=(14, 10))\n\n# Bar chart for IA categories\nax = axes[0, 0]\ncolors_ia = ['#2ECC71', '#3498DB', '#F39C12', '#E74C3C', '#8E44AD']\ncounts = [ia_category_counts.get(label, 0) for label in ia_labels]\nbars = ax.bar(range(len(ia_labels)), counts, color=colors_ia, alpha=0.7, edgecolor='black', linewidth=2)\nax.set_xticks(range(len(ia_labels)))\nax.set_xticklabels(ia_labels, rotation=15, ha='right')\nax.set_ylabel('Number of GO Terms', fontweight='bold')\nax.set_title('Distribution of GO Terms by IA Score Range', fontweight='bold')\n\n# Add labels\nfor bar, val, pct in zip(bars, counts, [c / len(ia_scores) * 100 for c in counts]):\n    ax.text(bar.get_x() + bar.get_width()/2, bar.get_height(),\n           f'{int(val):,}\\n({pct:.1f}%)', ha='center', va='bottom', fontweight='bold', fontsize=9)\n\n# Histogram of all IA scores\nax = axes[0, 1]\nax.hist(ia_scores['IA_score'].values, bins=50, color='indianred', alpha=0.7, edgecolor='black')\nax.axvline(ia_scores['IA_score'].mean(), color='red', linestyle='--', linewidth=2, label=f'Mean: {ia_scores[\"IA_score\"].mean():.2f}')\nax.axvline(ia_scores['IA_score'].median(), color='green', linestyle='--', linewidth=2, label=f'Median: {ia_scores[\"IA_score\"].median():.2f}')\nax.set_xlabel('IA Score', fontweight='bold')\nax.set_ylabel('Frequency', fontweight='bold')\nax.set_title('Histogram of All IA Scores', fontweight='bold')\nax.legend()\nax.set_yscale('log')\n\n# Pie chart\nax = axes[1, 0]\npie_labels = [f'{label}\\n({ia_type_labels[i]})' for i, label in enumerate(ia_labels)]\nax.pie(counts, labels=pie_labels, autopct='%1.1f%%', colors=colors_ia, startangle=90,\n      textprops={'fontsize': 9, 'fontweight': 'bold'})\nax.set_title('IA Score Categories Distribution', fontweight='bold')\n\n# Cumulative distribution\nax = axes[1, 1]\nsorted_ia = np.sort(ia_scores['IA_score'].values)\ncumulative_ia = np.arange(1, len(sorted_ia) + 1) / len(sorted_ia) * 100\nax.plot(sorted_ia, cumulative_ia, linewidth=2.5, color='darkblue')\nax.axvline(3, color='orange', linestyle='--', alpha=0.7, linewidth=2, label='Rare threshold (IA=3)')\nax.fill_between(sorted_ia, cumulative_ia, alpha=0.2, color='blue')\nax.set_xlabel('IA Score', fontweight='bold')\nax.set_ylabel('Cumulative % of GO Terms', fontweight='bold')\nax.set_title('Cumulative Distribution of IA Scores', fontweight='bold')\nax.legend()\nax.grid(True, alpha=0.3)\n\nplt.tight_layout()\nplt.show()\n\nprint(f\"\\n✓ IA score distribution visualized\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:25.851006Z","iopub.execute_input":"2025-12-17T10:29:25.851595Z","iopub.status.idle":"2025-12-17T10:29:27.142373Z","shell.execute_reply.started":"2025-12-17T10:29:25.851567Z","shell.execute_reply":"2025-12-17T10:29:27.141673Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 8. Species Distribution","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SECTION 6: PHÂN BỔ THEO LOÀI (SPECIES)\")\nprint(\"=\"*80)\n\n# Count proteins per species\nspecies_counts = train_taxonomy['TaxonomyID'].value_counts()\n\n# Map taxonomy to species names\ntaxonomy_map = dict(zip(test_taxonomy['ID'], test_taxonomy['Species']))\n\n# Get species names\nspecies_names = []\nfor tax_id in species_counts.index[:10]:  # Top 10\n    name = taxonomy_map.get(tax_id, f'TaxID:{tax_id}')\n    species_names.append(name)\n\n# Top 10 species\nprint(\"\\nTop 10 Species by Protein Count:\")\ntop_species_df = pd.DataFrame({\n    'Rank': range(1, 11),\n    'TaxonomyID': species_counts.index[:10].values,\n    'Species Name': species_names,\n    'Protein Count': species_counts.values[:10],\n    '%': (species_counts.values[:10] / species_counts.sum() * 100).round(1)\n})\n\nprint(top_species_df.to_string(index=False))\n\n# Other species summary\nother_count = species_counts[10:].sum()\nn_other_species = len(species_counts) - 10\nprint(f\"\\n... ({n_other_species} other species): {other_count:,} proteins ({other_count/species_counts.sum()*100:.1f}%)\")\nprint(f\"\\nTotal unique species: {len(species_counts)}\")\n\n# Visualize\nfig, axes = plt.subplots(1, 2, figsize=(16, 6))\n\n# Top 15 species bar chart\nax = axes[0]\ntop_15_taxa = species_counts.head(15)\ntop_15_names = [taxonomy_map.get(tid, f'TaxID:{tid}') for tid in top_15_taxa.index]\nbars = ax.barh(range(len(top_15_taxa)), top_15_taxa.values, color='teal', alpha=0.7, edgecolor='black')\nax.set_yticks(range(len(top_15_taxa)))\nax.set_yticklabels(top_15_names, fontsize=9)\nax.set_xlabel('Protein Count', fontweight='bold')\nax.set_title('Top 15 Species by Number of Proteins', fontweight='bold')\nax.invert_yaxis()\n\n# Add value labels\nfor i, (bar, val) in enumerate(zip(bars, top_15_taxa.values)):\n    pct = val / species_counts.sum() * 100\n    ax.text(val, bar.get_y() + bar.get_height()/2, f' {val:,} ({pct:.1f}%)',\n           ha='left', va='center', fontweight='bold', fontsize=9)\n\n# Pie chart for top 10\nax = axes[1]\ntop_10_taxa = species_counts.head(10)\ntop_10_names = [taxonomy_map.get(tid, f'TaxID:{tid}').split()[0] for tid in top_10_taxa.index]  # First word\nothers = species_counts[10:].sum()\n\npie_values = list(top_10_taxa.values) + [others]\npie_labels = top_10_names + ['Others']\n\ncolors_species = plt.cm.Set3(np.linspace(0, 1, len(pie_values)))\nwedges, texts, autotexts = ax.pie(pie_values, labels=pie_labels, autopct='%1.1f%%',\n                                    colors=colors_species, startangle=90,\n                                    textprops={'fontsize': 9, 'fontweight': 'bold'})\nax.set_title('Species Distribution (Top 10 + Others)', fontweight='bold')\n\nplt.tight_layout()\nplt.show()\n\nprint(f\"\\n✓ Species distribution visualized\")","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-12-17T10:29:27.143117Z","iopub.execute_input":"2025-12-17T10:29:27.1434Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"## 9. Summary Report","metadata":{}},{"cell_type":"code","source":"print(\"\\n\" + \"=\"*80)\nprint(\"SUMMARY REPORT - CAFA06 DATA OVERVIEW\")\nprint(\"=\"*80)\n\nprint(f\"\"\"\n╔═══════════════════════════════════════════════════════════════════════════════╗\n║                         DATA STRUCTURE SUMMARY                                ║\n╚═══════════════════════════════════════════════════════════════════════════════╝\n\n📊 DATASET SIZE:\n   • Training proteins:           {n_train_proteins:>10,}\n   • Test proteins:               {n_test_proteins:>10,}\n   • Total annotations:           {n_annotations:>10,}\n   • Unique GO terms:             {n_unique_go_terms:>10,}\n   • Unique species:              {n_unique_species:>10}\n\n🎯 ASPECT BREAKDOWN (GO Annotation Types):\n   • Cellular Component (C):      {aspect_counts['C']:>10,} ({aspect_counts['C']/aspect_counts.sum()*100:>5.1f}%)\n   • Molecular Function (F):      {aspect_counts['F']:>10,} ({aspect_counts['F']/aspect_counts.sum()*100:>5.1f}%)\n   • Biological Process (P):      {aspect_counts['P']:>10,} ({aspect_counts['P']/aspect_counts.sum()*100:>5.1f}%)\n\n📈 GO TERM STATISTICS:\n   • Most common term:            {go_counts.index[0]} ({go_counts.iloc[0]:,} proteins)\n   • Mean annotations/term:       {go_counts.mean():>10.1f}\n   • Median annotations/term:     {go_counts.median():>10.0f}\n   • Std Dev:                     {go_counts.std():>10.1f}\n\n🧬 PROTEIN ANNOTATION STATISTICS:\n   • Mean annotations/protein:    {protein_counts.mean():>10.1f}\n   • Median annotations/protein:  {protein_counts.median():>10.0f}\n   • Max annotations/protein:     {protein_counts.max():>10}\n   • Std Dev:                     {protein_counts.std():>10.1f}\n\n🔍 IA SCORE ANALYSIS (Term Rarity):\n   • Very Common (IA 0.0-1.0):   {ia_category_counts[0]:>10,} ({ia_category_counts[0]/len(ia_scores)*100:>5.1f}%)\n   • Common (IA 1.0-2.0):         {ia_category_counts[1]:>10,} ({ia_category_counts[1]/len(ia_scores)*100:>5.1f}%)\n   • Semi-Rare (IA 2.0-3.0):      {ia_category_counts[2]:>10,} ({ia_category_counts[2]/len(ia_scores)*100:>5.1f}%)\n   • Rare (IA 3.0-4.0):           {ia_category_counts[3]:>10,} ({ia_category_counts[3]/len(ia_scores)*100:>5.1f}%)\n   • Very Rare (IA 4.0-6.0):      {ia_category_counts[4]:>10,} ({ia_category_counts[4]/len(ia_scores)*100:>5.1f}%)\n   • Mean IA score:               {ia_scores['IA_score'].mean():>10.2f}\n   • Median IA score:             {ia_scores['IA_score'].median():>10.2f}\n\n🌍 TOP 5 SPECIES:\n   1. {taxonomy_map.get(species_counts.index[0], 'Unknown'):<40} {species_counts.iloc[0]:>7,} proteins\n   2. {taxonomy_map.get(species_counts.index[1], 'Unknown'):<40} {species_counts.iloc[1]:>7,} proteins\n   3. {taxonomy_map.get(species_counts.index[2], 'Unknown'):<40} {species_counts.iloc[2]:>7,} proteins\n   4. {taxonomy_map.get(species_counts.index[3], 'Unknown'):<40} {species_counts.iloc[3]:>7,} proteins\n   5. {taxonomy_map.get(species_counts.index[4], 'Unknown'):<40} {species_counts.iloc[4]:>7,} proteins\n\n\"\"\")\n\nprint(\"=\"*80)\nprint(\"✓ Data analysis completed successfully!\")\nprint(\"=\"*80)","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}