{"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":"# Outcomes: \n\n## Samples counts:\n\nTrain and test about 140 000 sequences each.\n\n## Labels\n\nIt is multi-label, multi-class task.\nThat means one sample should be assigned MANY labels. \nMoreover there are certain weights for the labels - which can be zero - meaning that prediction of those - will not affect evaluation. (We will analyse that part later). \n\n    unique terms(labels) in the train: 31466\n    unique terms(labels) by the Gene Ontology aspect: \n    BPO    21285\n    CCO     2957\n    MFO     7224\n\nIn average one label is assigned to 170 samples, but median is only 8.\nAnd moreover: \n\n    12 terms(labels) occured more than in  50000 samples (train)\n    82 terms(labels) occured more than in  10000 samples (train)\n    163 terms(labels) occured more than in  5000 samples (train)\n    771 terms(labels) occured more than in  1000 samples (train)\n\n\nIn train 37 labels are assigned in average to each sample (median 24).\n\n    That 37 can be splitted into 3 GO-subcategories:\n    BPO    24.589317 - Biological Process - aspect of Gene Ontology \n    CCO     8.408089 - Cell Component - aspect of  Gene Ontology\n    MFO     4.710951 - Molecular function - aspect of Gene Ontology \n\nPrevious winner had the most advantage in MFO - may be focus on that? but it is the smallest category. \nSo not clear what to suggest. \nhttps://www.kaggle.com/competitions/cafa-5-protein-function-prediction/discussion/402645#2227469\n\n    37 is average\n    Here we look on the distribution of how many lables is assigned to each sample in the train: \n    9875 count samples with more than  100   terms(=labels) in the train\n    1553 count samples with more than  200   terms(=labels) in the train\n    25 count samples with more than  500   terms(=labels) in the train\n\n    12972 count samples with less than  5   terms(=labels) in the train\n    32939 count samples with less than  10   terms(=labels) in the train\n    105782 count samples with less than  50   terms(=labels) in the train\n\n\n## Taxons:\n\nThere are 90 organisms (taxons) in the TEST and more than 3000 in the TRAIN.\n\nNevertheless top represented taxons well agree in the train and test:\n\n    Count intersection of top N taxons in train and test:\n    N= 10 n_common: 8\n    N= 25 n_common: 20\n    N= 50 n_common: 41\n    N= 90 n_common: 70\n\nNOT(!) all taxons from the TEST presented in the TRAIN  - only 71. \n\n    Count taxons with more than 1000 sample in train: 18\n    Count taxons with more than 100 sample in train: 48\n\n    Top 5 presented taxons (=ogranisms = species ) in the TRAIN:\n    taxonomyID\t\tCount  Train\tSpecies\n    9606\t\t\t25125\t\t\thomo sapiens[All Names]\n    3702\t\t\t14461\t\t\tArabidopsis thaliana[All Names]  --  a small plant from the mustard family \n    10090\t\t\t14384\t\t\tmus musculus[All Names]\n    7955\t\t\t12671\t\t\tDanio rerio[All Names]\n    7227\t\t\t12020\t\t\tDrosophila melanogaster[All Names]\n\n\n\n","metadata":{}},{"cell_type":"markdown","source":"# Draft code below. To be updated\n\nSome borrowed from the previous EDA: \n\nhttps://www.kaggle.com/code/leonidkulyk/eda-cafa5-pfp-interactive-dags-plotly\n\nhttps://www.kaggle.com/code/mpwolke/cafa-5-protein-prediction\n\nhttps://www.kaggle.com/code/thedrcat/cafa-eda\n\nThanks to the authors ! \n","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 time\nt0start = time.time()\n\nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\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-07-21T09:25:31.700059Z","iopub.execute_input":"2023-07-21T09:25:31.700821Z","iopub.status.idle":"2023-07-21T09:25:33.092721Z","shell.execute_reply.started":"2023-07-21T09:25:31.700762Z","shell.execute_reply":"2023-07-21T09:25:33.091222Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!pip install obonet\n!pip install pyvis\n\nimport os\nimport json\nfrom PIL import Image\nfrom typing import Dict\nfrom collections import Counter\n\nimport random\nimport cv2\nimport obonet\nimport networkx\nimport pandas as pd\nimport numpy as np\nimport plotly.express as px\nimport plotly.graph_objects as go\nimport matplotlib.pyplot as plt\nimport matplotlib.patches as mpatch\nimport seaborn  as sns\nfrom Bio import SeqIO\nfrom pyvis.network import Network","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:25:33.095886Z","iopub.execute_input":"2023-07-21T09:25:33.097015Z","iopub.status.idle":"2023-07-21T09:26:07.735078Z","shell.execute_reply.started":"2023-07-21T09:25:33.096958Z","shell.execute_reply":"2023-07-21T09:26:07.73382Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# train_taxonomy.tsv","metadata":{}},{"cell_type":"code","source":"fn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_taxonomy.tsv'\ndf = pd.read_csv(fn, sep = '\\t')\ndisplay(df)\nv_taxonomyID_value_counts_train = df['taxonomyID'].value_counts()  \nv_taxonomyID_value_counts_train.name = 'Count Taxons in Train'\ndisplay(v_taxonomyID_value_counts_train.head(30) )\nv_taxonomyID_value_counts_train.nunique()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:07.736899Z","iopub.execute_input":"2023-07-21T09:26:07.737395Z","iopub.status.idle":"2023-07-21T09:26:07.938132Z","shell.execute_reply.started":"2023-07-21T09:26:07.73734Z","shell.execute_reply":"2023-07-21T09:26:07.936501Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m = v_taxonomyID_value_counts_train > 100\nprint('Count taxons with more than 100 sample in train:',  m.sum() )\nm = v_taxonomyID_value_counts_train > 1000\nprint('Count taxons with more than 1000 sample in train:',  m.sum() )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:07.94322Z","iopub.execute_input":"2023-07-21T09:26:07.943704Z","iopub.status.idle":"2023-07-21T09:26:07.953347Z","shell.execute_reply.started":"2023-07-21T09:26:07.943644Z","shell.execute_reply":"2023-07-21T09:26:07.951902Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# testsuperset-taxon-list.tsv","metadata":{}},{"cell_type":"code","source":"from pathlib import Path\npath = Path('../input/cafa-5-protein-function-prediction')\n!head {path}/'Test (Targets)/testsuperset-taxon-list.tsv'","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:07.955221Z","iopub.execute_input":"2023-07-21T09:26:07.955868Z","iopub.status.idle":"2023-07-21T09:26:09.086264Z","shell.execute_reply.started":"2023-07-21T09:26:07.955823Z","shell.execute_reply":"2023-07-21T09:26:09.084666Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Information on taxons:  Id, Ogranism info, count in the train ","metadata":{}},{"cell_type":"code","source":"fn = '../input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset-taxon-list.tsv'\ndf = pd.read_csv(fn, sep='\\t', error_bad_lines=False, encoding= 'unicode_escape')\ndisplay(df)\n# display(df['Species'].value_counts().head(30) )\nprint(df['Species'].nunique() )\ns =   set( v_taxonomyID_value_counts_train.index ) & set(df ['ID'] )\nprint('Intersection of taxons on train and test - NOT complete :  ', len(s), 'out of', len(  set(df ['ID'] ) ) )\nprint(  s ) \n# s =   set( v_taxonomyID_value_counts_train.index ) | set(df ['ID'] )\n# print(len(s) ) \ndf_taxons_inf = v_taxonomyID_value_counts_train.to_frame().join( df.set_index('ID') , how = 'inner'  )#.head(50)\ndf_taxons_inf.index.name = 'taxonomyID'\nprint()\nprint('Common taxons for train and test:')\ndisplay( df_taxons_inf.head(50) )\ndisplay( df_taxons_inf.tail(21) )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:09.08935Z","iopub.execute_input":"2023-07-21T09:26:09.091001Z","iopub.status.idle":"2023-07-21T09:26:09.163863Z","shell.execute_reply.started":"2023-07-21T09:26:09.090937Z","shell.execute_reply":"2023-07-21T09:26:09.16221Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# testsuperset.fasta","metadata":{}},{"cell_type":"code","source":"!head {path}/'Test (Targets)/testsuperset.fasta'","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:09.166024Z","iopub.execute_input":"2023-07-21T09:26:09.16642Z","iopub.status.idle":"2023-07-21T09:26:10.272494Z","shell.execute_reply.started":"2023-07-21T09:26:09.166382Z","shell.execute_reply":"2023-07-21T09:26:10.270751Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nprint(\"Sequence example:\\n\\n\", next(iter(SeqIO.parse(fn, \"fasta\"))))\nsequences = SeqIO.parse(fn, \"fasta\")\nnum_sequences = sum(1 for seq in sequences)\nprint()\nprint(\"Number of sequences in super-test:\", num_sequences)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:10.274773Z","iopub.execute_input":"2023-07-21T09:26:10.275194Z","iopub.status.idle":"2023-07-21T09:26:12.952258Z","shell.execute_reply.started":"2023-07-21T09:26:10.275147Z","shell.execute_reply":"2023-07-21T09:26:12.950923Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from Bio import SeqIO\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\n\nseq = next(iter(sequences))\ngb = seq\nprint('\\nLength of Sequence:')\nprint(len(gb.seq))\n\nprint('\\nRecord ID:')\nprint(gb.id)\n\nprint('\\nName:')\nprint(gb.name)\n\nprint('\\nDescription:')\nprint(gb.description)\n\n# Annotations \nprint('\\nNumber of Annotations:')\nprint(len(gb.annotations))\n\n# Features \nprint('\\nNumber of Features:')\nprint(len(gb.features))","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:12.954333Z","iopub.execute_input":"2023-07-21T09:26:12.954902Z","iopub.status.idle":"2023-07-21T09:26:12.96612Z","shell.execute_reply.started":"2023-07-21T09:26:12.954861Z","shell.execute_reply":"2023-07-21T09:26:12.964892Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fn = '../input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset-taxon-list.tsv'\ndf = pd.read_csv(fn, sep='\\t', error_bad_lines=False, encoding= 'unicode_escape', index_col = 0)\ndisplay(df)\ndict_num2taxon_name = dict( df['Species'] )\nstr(dict_num2taxon_name)[:100]","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:12.972324Z","iopub.execute_input":"2023-07-21T09:26:12.972823Z","iopub.status.idle":"2023-07-21T09:26:13.014079Z","shell.execute_reply.started":"2023-07-21T09:26:12.97278Z","shell.execute_reply":"2023-07-21T09:26:13.012966Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfrom Bio import SeqIO\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\n\ndf_test_seq = pd.DataFrame(); i = 0;\n\nsequences = SeqIO.parse(fn, \"fasta\")\nl = [str(seq.description).split('\\t')[0] for seq in sequences]\ndf_test_seq['Id'] = l\n\nsequences = SeqIO.parse(fn, \"fasta\")\nl = [str(seq.description).split('\\t')[-1] for seq in sequences]\ndf_test_seq['Taxon'] = l\n\ndict_num2taxon_name[8667] = 'Australian taipan'\nsequences = SeqIO.parse(fn, \"fasta\")\nl = [dict_num2taxon_name[ float(str(seq.description).split('\\t')[-1]) ] for seq in sequences ]\ndf_test_seq['Taxon Name'] = l\n\nsequences = SeqIO.parse(fn, \"fasta\")\nl = [len(str(seq.seq)) for seq in sequences ]\ndf_test_seq['Sequence Len'] = l\n\nsequences = SeqIO.parse(fn, \"fasta\")\nl = [str(seq.seq) for seq in sequences ]\ndf_test_seq['Sequence'] = l\ndf_test_seq    ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:13.016035Z","iopub.execute_input":"2023-07-21T09:26:13.01734Z","iopub.status.idle":"2023-07-21T09:26:26.40598Z","shell.execute_reply.started":"2023-07-21T09:26:13.017295Z","shell.execute_reply":"2023-07-21T09:26:26.404595Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df_test_seq.tail(50)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:26.407966Z","iopub.execute_input":"2023-07-21T09:26:26.408813Z","iopub.status.idle":"2023-07-21T09:26:26.434721Z","shell.execute_reply.started":"2023-07-21T09:26:26.40876Z","shell.execute_reply":"2023-07-21T09:26:26.433378Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\ndf_test_seq.to_csv('df_test_seq_info.csv')","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:26.436535Z","iopub.execute_input":"2023-07-21T09:26:26.437045Z","iopub.status.idle":"2023-07-21T09:26:30.202861Z","shell.execute_reply.started":"2023-07-21T09:26:26.436996Z","shell.execute_reply":"2023-07-21T09:26:30.20166Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Example to work with Seq object","metadata":{}},{"cell_type":"code","source":"fn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\n\nsequences = SeqIO.parse(fn, \"fasta\")\nl = []\nseq = next(iter(sequences))\n\nprint('seq example:\\n', seq); print()\nl = [str(t) for t in dir(seq) if not str(t).startswith('_')]\nprint('seq attributes:', l);print()\nfor a in l:\n    print(a)\n    print(getattr(seq,a))\n    \n    \n# dir(seq)    \n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:30.20422Z","iopub.execute_input":"2023-07-21T09:26:30.204568Z","iopub.status.idle":"2023-07-21T09:26:30.219279Z","shell.execute_reply.started":"2023-07-21T09:26:30.204534Z","shell.execute_reply":"2023-07-21T09:26:30.21773Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"list_taxons_test = []\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\nfor seq in sequences:\n    t = str(getattr(seq,'description')).split('\\t')[-1]\n    list_taxons_test.append(t)\n\ndf_taxons_info_test = pd.Series(list_taxons_test).value_counts().to_frame()\ndf_taxons_info_test.columns = ['Count in Test']\ndf_taxons_info_test.index = [ int(t) for t in df_taxons_info_test.index ]\ndf_taxons_info_test.index.name = 'ID'\nfn = '../input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset-taxon-list.tsv'\ndf = pd.read_csv(fn, sep='\\t', error_bad_lines=False, encoding= 'unicode_escape').set_index('ID')\n#display(df)\ndf_taxons_info_test =  df_taxons_info_test.join( df)\ndisplay(df_taxons_info_test.head(50))\ndisplay(df_taxons_info_test.tail(40))\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:30.221048Z","iopub.execute_input":"2023-07-21T09:26:30.221432Z","iopub.status.idle":"2023-07-21T09:26:32.981697Z","shell.execute_reply.started":"2023-07-21T09:26:30.221395Z","shell.execute_reply":"2023-07-21T09:26:32.980375Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimport numpy as np\nimport pandas as pd\nfrom Bio import SeqIO\n\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\nl = [len(seq) for seq in sequences ] \nprint(pd.Series(l).describe() )\n\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nplt.figure(figsize = (20,3) )\nplt.hist(l, bins = 100 )\nplt.title('test lengths', fontsize = 20 )\nplt.show()\nplt.figure(figsize = (20,3) )\nplt.hist([t for t in l if t < 3000], bins = 100 )\nplt.title('test lengths', fontsize = 20 )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:32.983487Z","iopub.execute_input":"2023-07-21T09:26:32.984422Z","iopub.status.idle":"2023-07-21T09:26:38.208381Z","shell.execute_reply.started":"2023-07-21T09:26:32.984367Z","shell.execute_reply":"2023-07-21T09:26:38.206914Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('Count intersection of top N taxons in train and test:')\nfor N in [10,25,50,90]:\n    s = set(df_taxons_info_test.index[:N] ) & set( df_taxons_inf.index[:N]  )\n    print('N=', N, 'n_common:',  len(s))\n\nprint()\nprint('Top taxons(organisms) - left - test, right - train:')\nd = pd.concat( [df_taxons_info_test.reset_index().head(50) , df_taxons_inf.reset_index().head(50) ], axis = 1  )\nd","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:38.210315Z","iopub.execute_input":"2023-07-21T09:26:38.210842Z","iopub.status.idle":"2023-07-21T09:26:38.249284Z","shell.execute_reply.started":"2023-07-21T09:26:38.210786Z","shell.execute_reply":"2023-07-21T09:26:38.247611Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimport numpy as np\nimport pandas as pd\nfrom Bio import SeqIO\n\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\n# fn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\nl = np.array( [len(seq) for seq in sequences ] ) \nsequences = SeqIO.parse(fn, \"fasta\")\nl2 = np.array( [seq.id for seq in sequences ] ) \nprint(pd.Series(l).describe(percentiles = [0.001, 0.05, 0.01,0.1, 0.25, 0.5,0.75, 0.9, 0.95,  0.99,0.999]) )\nfor k in [1000, 2000, 3000, 5000, 10_000,15_000, 20_000,30_000]:\n    print('Len GE:',k, 'Count:',  (l>k).sum() )\nset_prev = set()\nfor k in [30_000,20_000,10_000, 8000 ]:\n    print('Len GE:',k ,set(l2[l>k])-set_prev )\n    set_prev = set(l2[l>k])\n    \nimport matplotlib.pyplot as plt\nimport seaborn as sns\nplt.figure(figsize = (20,3) )\nplt.hist(l, bins = 100 )\nplt.title('test lengths', fontsize = 20 )\nplt.show()\nplt.figure(figsize = (20,3) )\nplt.hist([t for t in l if t < 3000], bins = 100 )\nplt.title('test lengths', fontsize = 20 )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:38.251218Z","iopub.execute_input":"2023-07-21T09:26:38.252213Z","iopub.status.idle":"2023-07-21T09:26:45.173897Z","shell.execute_reply.started":"2023-07-21T09:26:38.252157Z","shell.execute_reply":"2023-07-21T09:26:45.172202Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# train_sequences","metadata":{}},{"cell_type":"code","source":"fn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nprint(\"Sequence example:\\n\\n\", next(iter(SeqIO.parse(fn, \"fasta\"))))\nsequences = SeqIO.parse(fn, \"fasta\")\nnum_sequences = sum(1 for seq in sequences)\nprint()\nprint(\"Number of sequences in super-test:\", num_sequences)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:45.176285Z","iopub.execute_input":"2023-07-21T09:26:45.176844Z","iopub.status.idle":"2023-07-21T09:26:48.032867Z","shell.execute_reply.started":"2023-07-21T09:26:45.176787Z","shell.execute_reply":"2023-07-21T09:26:48.031328Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimport numpy as np \nimport pandas as pd \nfrom Bio import SeqIO\n\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\n# l =  [set(list(seq.seq)) for seq in sequences ] \ns = set()\nfor seq in sequences:\n    s = set(list(seq.seq)) | s\nprint(len(s), s )\ns1 = s\n\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\n# l =  [set(list(seq.seq)) for seq in sequences ] \ns = set()\nfor seq in sequences:\n    s = set(list(seq.seq)) | s\nprint(len(s), s )\ns2 = s \nprint(s1-s2, s2-s1 )\n\nletters = ['Y', 'K', 'Z', 'U', 'R', 'N', 'D', 'B', 'S', 'P', 'W', 'T', 'G', 'H', 'O', 'L', 'E', 'C', 'A', 'I', 'X', 'F', 'V', 'M', 'Q']\n\ndict_freq_train = {}\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nfor k in letters:\n    sequences = SeqIO.parse(fn, \"fasta\")\n    l =  [str(seq.seq).count(k) for seq in sequences ] \n    dict_freq_train[k] = np.sum(l)\ndict_freq_test = {}\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nfor k in letters:\n    sequences = SeqIO.parse(fn, \"fasta\")\n    l =  [str(seq.seq).count(k) for seq in sequences ] \n    dict_freq_test[k] = np.sum(l)\n    \nd = pd.concat([ pd.Series(dict_freq_train).to_frame().reset_index() , pd.Series(dict_freq_test).to_frame().reset_index() ], axis = 1) \nd.columns = ['Train Amino Acid', 'Train Count', 'Test Amino Acid', 'Test Count', ]\ndisplay(d.sort_values('Train Count'))","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:26:48.034787Z","iopub.execute_input":"2023-07-21T09:26:48.035402Z","iopub.status.idle":"2023-07-21T09:29:23.40303Z","shell.execute_reply.started":"2023-07-21T09:26:48.035363Z","shell.execute_reply":"2023-07-21T09:29:23.40156Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import pandas as pd\n\ndata = {'Amino Acid': ['Alanine', 'Arginine', 'Asparagine', 'Aspartic acid', 'Cysteine', 'Glutamine', 'Glutamic acid', 'Glycine', 'Histidine', 'Isoleucine', 'Leucine', 'Lysine', 'Methionine', 'Phenylalanine', 'Proline', 'Pyrrolysine', 'Serine', 'Selenocysteine', 'Threonine', 'Tryptophan', 'Tyrosine', 'Valine', 'Aspartic acid or Asparagine', 'Glutamic acid or Glutamine', 'Any amino acid', 'Leucine or Isoleucine'],\n        '1-letter Code': ['A', 'R', 'N', 'D', 'C', 'Q', 'E', 'G', 'H', 'I', 'L', 'K', 'M', 'F', 'P', 'O', 'S', 'U', 'T', 'W', 'Y', 'V', 'B', 'Z', 'X', 'J'],\n        '3-letter Code': ['Ala', 'Arg', 'Asn', 'Asp', 'Cys', 'Gln', 'Glu', 'Gly', 'His', 'Ile', 'Leu', 'Lys', 'Met', 'Phe', 'Pro', 'Pyl', 'Ser', 'Sec', 'Thr', 'Trp', 'Tyr', 'Val', 'Asx', 'Glx', 'Xaa', 'Xle']}\n        \ndf = pd.DataFrame(data)\nd2 = pd.merge( d, df, how = 'left', left_on = 'Test Amino Acid' , right_on = '1-letter Code' )\nd2 = d2.sort_values('Test Count')\nd2","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:23.405032Z","iopub.execute_input":"2023-07-21T09:29:23.405537Z","iopub.status.idle":"2023-07-21T09:29:23.443498Z","shell.execute_reply.started":"2023-07-21T09:29:23.405485Z","shell.execute_reply":"2023-07-21T09:29:23.442386Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d2.to_csv('amino_acid_counts.csv')","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:23.445122Z","iopub.execute_input":"2023-07-21T09:29:23.445512Z","iopub.status.idle":"2023-07-21T09:29:23.46415Z","shell.execute_reply.started":"2023-07-21T09:29:23.445473Z","shell.execute_reply":"2023-07-21T09:29:23.462702Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d2[['Train Count', 'Test Count']].diff()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:23.469124Z","iopub.execute_input":"2023-07-21T09:29:23.469626Z","iopub.status.idle":"2023-07-21T09:29:23.498159Z","shell.execute_reply.started":"2023-07-21T09:29:23.46958Z","shell.execute_reply":"2023-07-21T09:29:23.496655Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimport numpy as np\nimport pandas as pd\nfrom Bio import SeqIO\n\n# fn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\nl = np.array( [len(seq) for seq in sequences ] ) \nsequences = SeqIO.parse(fn, \"fasta\")\nl2 = np.array( [seq.id for seq in sequences ] ) \nprint(pd.Series(l).describe(percentiles = [0.001, 0.05, 0.01,0.1, 0.25, 0.5,0.75, 0.9, 0.95,  0.99,0.999]) )\nfor k in [1000, 2000, 3000, 5000, 10_000,15_000, 20_000,30_000]:\n    print('Len GE:', k, 'Count:',  (l>k).sum() )\nset_prev = set()\nfor k in [30_000,20_000,10_000, 8000 ]:\n    print('Len GE:', k, set(l2[l>k])-set_prev )\n    set_prev = set(l2[l>k])\n    \nimport matplotlib.pyplot as plt\nimport seaborn as sns\nplt.figure(figsize = (20,3) )\nplt.hist(l, bins = 100 )\nplt.title('train lengths', fontsize = 20 )\nplt.show()\nplt.figure(figsize = (20,3) )\nplt.hist([t for t in l if t < 3000], bins = 100 )\nplt.title('train lengths', fontsize = 20 )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:23.500054Z","iopub.execute_input":"2023-07-21T09:29:23.500456Z","iopub.status.idle":"2023-07-21T09:29:30.754964Z","shell.execute_reply.started":"2023-07-21T09:29:23.500416Z","shell.execute_reply":"2023-07-21T09:29:30.753616Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfrom Bio import SeqIO\nimport numpy as np\n\n# fn = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\n\nsequences = SeqIO.parse(fn, \"fasta\")\n\ni = 0; l = []\nll = []; ll2 = []\nfor seq in sequences: #   in range(5000):\n    ll.append( len(seq.seq) )\n    ll2.append(seq.id)\nll = np.array( ll )\nll2 = np.array(ll2)\nprint( len(ll) )\nv = 0\nfor i in range(20):\n    m = ll == i\n    v += m.sum()\n    print(i, m.sum(), v, ll2[m] )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:30.756537Z","iopub.execute_input":"2023-07-21T09:29:30.756964Z","iopub.status.idle":"2023-07-21T09:29:33.637831Z","shell.execute_reply.started":"2023-07-21T09:29:30.756926Z","shell.execute_reply":"2023-07-21T09:29:33.636099Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfrom Bio import SeqIO\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\n\nseq = next(iter(sequences))\ngb = seq\nprint('\\nLength of Sequence:')\nprint(len(gb.seq))\n\nprint('\\nRecord ID:')\nprint(gb.id)\n\nprint('\\nName:')\nprint(gb.name)\n\nprint('\\nDescription:')\nprint(gb.description)\n\n# Annotations \nprint('\\nNumber of Annotations:')\nprint(len(gb.annotations))\n\n# Features \nprint('\\nNumber of Features:')\nprint(len(gb.features))","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:33.639829Z","iopub.execute_input":"2023-07-21T09:29:33.640225Z","iopub.status.idle":"2023-07-21T09:29:33.651523Z","shell.execute_reply.started":"2023-07-21T09:29:33.640187Z","shell.execute_reply":"2023-07-21T09:29:33.649668Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# train_terms","metadata":{}},{"cell_type":"code","source":"!head {path}/'Train/train_terms.tsv'","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:33.653566Z","iopub.execute_input":"2023-07-21T09:29:33.654146Z","iopub.status.idle":"2023-07-21T09:29:34.777384Z","shell.execute_reply.started":"2023-07-21T09:29:33.654102Z","shell.execute_reply":"2023-07-21T09:29:34.776061Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fn =  str(path) + '/Train/train_terms.tsv'\ndf = pd.read_csv(fn , sep = '\\t')\ndf","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:34.791681Z","iopub.execute_input":"2023-07-21T09:29:34.792156Z","iopub.status.idle":"2023-07-21T09:29:38.76332Z","shell.execute_reply.started":"2023-07-21T09:29:34.792112Z","shell.execute_reply":"2023-07-21T09:29:38.762104Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('unique terms in the train:', df['term'].nunique() )\nprint('unique terms by the aspect: ' )\ndf.groupby('aspect')['term'].nunique()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:38.765051Z","iopub.execute_input":"2023-07-21T09:29:38.765439Z","iopub.status.idle":"2023-07-21T09:29:41.544295Z","shell.execute_reply.started":"2023-07-21T09:29:38.765382Z","shell.execute_reply":"2023-07-21T09:29:41.543165Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"len(df) / num_sequences","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:41.547078Z","iopub.execute_input":"2023-07-21T09:29:41.54758Z","iopub.status.idle":"2023-07-21T09:29:41.558572Z","shell.execute_reply.started":"2023-07-21T09:29:41.547528Z","shell.execute_reply":"2023-07-21T09:29:41.557222Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df['aspect'].value_counts()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:41.560616Z","iopub.execute_input":"2023-07-21T09:29:41.561145Z","iopub.status.idle":"2023-07-21T09:29:42.505154Z","shell.execute_reply.started":"2023-07-21T09:29:41.561091Z","shell.execute_reply":"2023-07-21T09:29:42.503643Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df['aspect'].value_counts()/ num_sequences","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:42.507239Z","iopub.execute_input":"2023-07-21T09:29:42.508168Z","iopub.status.idle":"2023-07-21T09:29:43.448716Z","shell.execute_reply.started":"2023-07-21T09:29:42.508123Z","shell.execute_reply":"2023-07-21T09:29:43.447401Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"v = df['term'].value_counts()\nprint('Statistics on terms(labels):')\nprint(v.describe() )\nfor k in [300,500,1000,1500,2000,3000, 5000,10000,]:\n    print(k, v.values[k],  v.index[k])\nplt.figure(figsize = (20,3))\nplt.plot(v.head(300).values)\nplt.xlabel('Rank of label (sorted by count of samples with such label)', fontsize = 20 )\nplt.ylabel('Count samples with such label', fontsize = 20  )\nplt.show()\n\nplt.figure(figsize = (20,3))\nplt.plot(range(300,3000),v.values[300:3000])\nplt.xlabel('Rank of label (sorted by count of samples with such label)', fontsize = 20 )\nplt.ylabel('Count samples with such label', fontsize = 20  )\nplt.show()\n\nplt.figure(figsize = (20,7))\nplt.title('Log10X Log10Y plot', fontsize = 20 )\nplt.plot(np.log10(range(100,20000)),np.log10(v.values[100:20000]) )\nplt.xlabel('Rank of label  (sorted by count of samples with such label) ', fontsize = 20 )\nplt.ylabel('Count samples with such label', fontsize = 20 )\nplt.show()\n\nplt.figure(figsize = (20,3))\nplt.hist(v[v<1000], bins = 200)\nplt.show()\n\nprint('Frequent terms(labels):')\nfor c in [50000, 10000, 5000, 1000]:\n    m = v > c\n    print(m.sum(), 'terms(labels) occured more than in ', c , 'samples (train)' )\n\nprint()\nprint('Rare terms(labels):')\nfor c in [2, 5, 10, 100]:\n    m = v < c\n    print(m.sum(), 'terms occured less than in ', c , 'samples (train)' )\n    \nprint()\nprint('Top frequent terms( labels ):')\ndisplay( v.head(50) )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:43.450263Z","iopub.execute_input":"2023-07-21T09:29:43.451137Z","iopub.status.idle":"2023-07-21T09:29:46.007379Z","shell.execute_reply.started":"2023-07-21T09:29:43.451094Z","shell.execute_reply":"2023-07-21T09:29:46.006023Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d = df.groupby('EntryID').count().sort_values('term', ascending = False)\nd ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:46.009147Z","iopub.execute_input":"2023-07-21T09:29:46.009551Z","iopub.status.idle":"2023-07-21T09:29:48.401172Z","shell.execute_reply.started":"2023-07-21T09:29:46.009502Z","shell.execute_reply":"2023-07-21T09:29:48.399778Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"d.describe()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:48.403042Z","iopub.execute_input":"2023-07-21T09:29:48.403416Z","iopub.status.idle":"2023-07-21T09:29:48.440112Z","shell.execute_reply.started":"2023-07-21T09:29:48.40338Z","shell.execute_reply":"2023-07-21T09:29:48.438747Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for c in [100,200,500]:\n    m = d['term'] > c\n    print(m.sum(), 'count samples with more than ', c, '  terms(=labels) in the train')\n    \nprint() ; print() ;    \nfor c in [5, 10,50]:\n    m = d['term'] < c\n    print(m.sum(), 'count samples with less than ', c, '  terms(=labels) in the train')    ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:48.442199Z","iopub.execute_input":"2023-07-21T09:29:48.442595Z","iopub.status.idle":"2023-07-21T09:29:48.457744Z","shell.execute_reply.started":"2023-07-21T09:29:48.442556Z","shell.execute_reply":"2023-07-21T09:29:48.456068Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import matplotlib.pyplot as plt\nimport seaborn as sns\nplt.figure(figsize = (20,3) )\nplt.plot( d['term'].values )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:48.459524Z","iopub.execute_input":"2023-07-21T09:29:48.460208Z","iopub.status.idle":"2023-07-21T09:29:48.720729Z","shell.execute_reply.started":"2023-07-21T09:29:48.460164Z","shell.execute_reply":"2023-07-21T09:29:48.719373Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.figure(figsize = (10,3) )\nplt.hist( d['term'].values, bins = 1000 )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:48.723494Z","iopub.execute_input":"2023-07-21T09:29:48.724798Z","iopub.status.idle":"2023-07-21T09:29:51.055276Z","shell.execute_reply.started":"2023-07-21T09:29:48.724752Z","shell.execute_reply":"2023-07-21T09:29:51.054162Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Lengths of sequences in the train - close look ","metadata":{}},{"cell_type":"code","source":"%%time \nfrom Bio import SeqIO\n\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'\nsequences = SeqIO.parse(fn, \"fasta\")\nl = []\nfor seq in sequences:\n    l.append(len(seq.seq))\ndisplay(pd.Series(l).describe() )\n\n%%time\nl = np.asarray(l)\nfor k in range(0,10):\n    m = l == k\n    print(k, m.sum() )\nsr = pd.Series(l).value_counts()\nsr = sr.sort_index()\ndisplay(sr )\n\nfig = plt.figure(figsize = (20,3) )\nplt.plot(sr.head(1500))\nplt.title('Train data',fontsize = 20 )\nplt.xlabel('Sequence length',fontsize = 20)\nplt.ylabel('Count',fontsize = 20)\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:51.056823Z","iopub.execute_input":"2023-07-21T09:29:51.057966Z","iopub.status.idle":"2023-07-21T09:29:53.811493Z","shell.execute_reply.started":"2023-07-21T09:29:51.057904Z","shell.execute_reply":"2023-07-21T09:29:53.810399Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# IA.txt - Information Accretion for each term. This is used to weight precision and recall (see Evaluation)","metadata":{}},{"cell_type":"code","source":"fn = '/kaggle/input/cafa-5-protein-function-prediction/IA.txt'\ndf_ia = pd.read_csv(fn , sep = '\\t', header = None)\ndisplay( df_ia )\n\ndisplay( df_ia.describe(percentiles = np.array([1,5,25,50,75, 95, 99 ])/100 ) )\n\nprint()\ns = (df_ia[1] == 0).sum()\nprint(s, 'count zero weight categories ', '%.1f - percent '%( s / len(df_ia)*100 ) )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:53.812848Z","iopub.execute_input":"2023-07-21T09:29:53.813708Z","iopub.status.idle":"2023-07-21T09:29:53.896037Z","shell.execute_reply.started":"2023-07-21T09:29:53.81367Z","shell.execute_reply":"2023-07-21T09:29:53.894683Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.hist(df_ia[1] , bins = 100  )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:53.89798Z","iopub.execute_input":"2023-07-21T09:29:53.899271Z","iopub.status.idle":"2023-07-21T09:29:54.318008Z","shell.execute_reply.started":"2023-07-21T09:29:53.899216Z","shell.execute_reply":"2023-07-21T09:29:54.316458Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"fn =  str(path) + '/Train/train_terms.tsv'\ndf_train_terms = pd.read_csv(fn , sep = '\\t')\ndf_train_terms\nm = (df_ia[1] == 0)\n\ns = set( df_train_terms['term'] ) & set( df_ia[m][0] )\nprint(len(s ), len(s )/len(set( df_train_terms['term'] )) ) #  ,   len(set( df_train_terms['term'] )  ) / len(df_ia),len(s) /  m.sum()   )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:54.319717Z","iopub.execute_input":"2023-07-21T09:29:54.320117Z","iopub.status.idle":"2023-07-21T09:29:59.025819Z","shell.execute_reply.started":"2023-07-21T09:29:54.320076Z","shell.execute_reply":"2023-07-21T09:29:59.024325Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m = df_ia[0].isin(  set( df_train_terms['term'] )  )\nplt.hist(df_ia[1][m] , bins = 100  )\nplt.show()\ndisplay( df_ia[m].describe(percentiles = np.array([1,5,25,50,75, 95, 99 ])/100 ) )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:29:59.02871Z","iopub.execute_input":"2023-07-21T09:29:59.030092Z","iopub.status.idle":"2023-07-21T09:30:00.409938Z","shell.execute_reply.started":"2023-07-21T09:29:59.030029Z","shell.execute_reply":"2023-07-21T09:30:00.408704Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nterms_frequencies = df_train_terms.groupby('term')['aspect'].count().sort_values(ascending = False)\nterms_frequencies","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:00.411687Z","iopub.execute_input":"2023-07-21T09:30:00.412372Z","iopub.status.idle":"2023-07-21T09:30:01.871042Z","shell.execute_reply.started":"2023-07-21T09:30:00.412329Z","shell.execute_reply":"2023-07-21T09:30:01.869496Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t = terms_frequencies.to_frame()\nt.columns = ['Frequency in Train']\nt2 = t.join(df_ia.set_index(0), how = 'left')\nt2.columns = [t2.columns[0], 'Weight Score']\nt2","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:01.875257Z","iopub.execute_input":"2023-07-21T09:30:01.875733Z","iopub.status.idle":"2023-07-21T09:30:01.918578Z","shell.execute_reply.started":"2023-07-21T09:30:01.875687Z","shell.execute_reply":"2023-07-21T09:30:01.917164Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"pd.concat([ t2.head(50).reset_index(), t2.tail(50).reset_index() ] , axis = 1 ) ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:01.920882Z","iopub.execute_input":"2023-07-21T09:30:01.921431Z","iopub.status.idle":"2023-07-21T09:30:01.954433Z","shell.execute_reply.started":"2023-07-21T09:30:01.921374Z","shell.execute_reply":"2023-07-21T09:30:01.952975Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\nplt.figure(figsize = (20,3))\nplt.plot( t2['Weight Score'].values[:1000] )\nplt.title(' IA - scoring weights for terms - terms ordered  (descendently)  by frequency in train ( the first 1000  )' , fontsize = 20  )\nplt.show()\nplt.figure(figsize = (20,3))\nplt.plot( t2['Weight Score'].values )\nplt.title('IA - scoring weights for terms - terms ordered (descendently) by frequency   in train', fontsize = 20 )\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:01.956319Z","iopub.execute_input":"2023-07-21T09:30:01.956917Z","iopub.status.idle":"2023-07-21T09:30:02.792552Z","shell.execute_reply.started":"2023-07-21T09:30:01.956788Z","shell.execute_reply":"2023-07-21T09:30:02.791003Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t2.corr(method = 'spearman')","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:02.794359Z","iopub.execute_input":"2023-07-21T09:30:02.794768Z","iopub.status.idle":"2023-07-21T09:30:02.818295Z","shell.execute_reply.started":"2023-07-21T09:30:02.79473Z","shell.execute_reply":"2023-07-21T09:30:02.816922Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t2.corr()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:02.819898Z","iopub.execute_input":"2023-07-21T09:30:02.820243Z","iopub.status.idle":"2023-07-21T09:30:02.834501Z","shell.execute_reply.started":"2023-07-21T09:30:02.82021Z","shell.execute_reply":"2023-07-21T09:30:02.832955Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimport seaborn as sns\nsns.scatterplot(x = t2.iloc[:,0],y = t2.iloc[:,1], )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:02.836154Z","iopub.execute_input":"2023-07-21T09:30:02.836648Z","iopub.status.idle":"2023-07-21T09:30:03.181747Z","shell.execute_reply.started":"2023-07-21T09:30:02.83659Z","shell.execute_reply":"2023-07-21T09:30:03.180319Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nimport seaborn as sns\nsns.scatterplot(x = t2.rank().iloc[:,0],y = t2.rank().iloc[:,1], )\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.183723Z","iopub.execute_input":"2023-07-21T09:30:03.184549Z","iopub.status.idle":"2023-07-21T09:30:03.559859Z","shell.execute_reply.started":"2023-07-21T09:30:03.184481Z","shell.execute_reply":"2023-07-21T09:30:03.558451Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"m = t2['Weight Score'] > 0\nt2[m].describe()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.561324Z","iopub.execute_input":"2023-07-21T09:30:03.561743Z","iopub.status.idle":"2023-07-21T09:30:03.591308Z","shell.execute_reply.started":"2023-07-21T09:30:03.561701Z","shell.execute_reply":"2023-07-21T09:30:03.589933Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for c1, c2 in [(1,0),(1,0.1),(1,1), (1,100),(5, 0),(5, 0.1), (5, 100),(10, 0),(10, 0.1), (10, 1),  \n               (100, 0),(100, 1), (100, 100), (1000, 0), (1000, 1), (1000, 2), (1000, 100),]: \n    m = (t2['Frequency in Train'] <= c1 ) & (t2['Weight Score'] <= c2)\n    print(c1,c2, m.sum() )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.592775Z","iopub.execute_input":"2023-07-21T09:30:03.593132Z","iopub.status.idle":"2023-07-21T09:30:03.616311Z","shell.execute_reply.started":"2023-07-21T09:30:03.593095Z","shell.execute_reply":"2023-07-21T09:30:03.61495Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"t2['CumSum Weight Score'] = t2['Weight Score'].cumsum()\nt2['Percent CumSum Weight Score'] = t2['CumSum Weight Score'] / t2['CumSum Weight Score'].iat[len(t2)-1] * 100 \nt2","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.617381Z","iopub.execute_input":"2023-07-21T09:30:03.617763Z","iopub.status.idle":"2023-07-21T09:30:03.64226Z","shell.execute_reply.started":"2023-07-21T09:30:03.617725Z","shell.execute_reply":"2023-07-21T09:30:03.641058Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"plt.plot(t2['Percent CumSum Weight Score'].values)\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.643816Z","iopub.execute_input":"2023-07-21T09:30:03.644166Z","iopub.status.idle":"2023-07-21T09:30:03.877557Z","shell.execute_reply.started":"2023-07-21T09:30:03.644132Z","shell.execute_reply":"2023-07-21T09:30:03.875916Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for c in [99, 95, 90 , 75]:\n    m  = t2['Percent CumSum Weight Score'] < c\n    print(c, m.sum() )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.879017Z","iopub.execute_input":"2023-07-21T09:30:03.879395Z","iopub.status.idle":"2023-07-21T09:30:03.890903Z","shell.execute_reply.started":"2023-07-21T09:30:03.879358Z","shell.execute_reply":"2023-07-21T09:30:03.889413Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Train terms - 2 ","metadata":{}},{"cell_type":"code","source":"fn =  str(path) + '/Train/train_terms.tsv'\ndf = pd.read_csv(fn , sep = '\\t')\ndf","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:03.892893Z","iopub.execute_input":"2023-07-21T09:30:03.893377Z","iopub.status.idle":"2023-07-21T09:30:06.859107Z","shell.execute_reply.started":"2023-07-21T09:30:03.893325Z","shell.execute_reply":"2023-07-21T09:30:06.857515Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time \ndisplay( df['term'].value_counts().head(20) )\ndisplay( df['term'].value_counts().tail(20) )\n\nfor t in range(1,11): # [1,2,3,4,5,10]:\n    print(t , (df['term'].value_counts() == t).sum() )\nv = df['term'].value_counts().value_counts()\nv = v.sort_index()\nplt.figure(figsize = (20,4))\nplt.plot(v.head(20),'*-')\nplt.title('count terms with give number of labels - show terms with least possible number of labels - 1,2,3...', fontsize = 20 )\nplt.xlabel('number of labels', fontsize = 20 )\nplt.grid()\nplt.show()\n\n# subontology_roots = {'BPO':'GO:0008150', #  Top2 frequent \n#                      'CCO':'GO:0005575', # Top frequent \n#                      'MFO':'GO:0003674'} # Top4 frequent","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:06.860777Z","iopub.execute_input":"2023-07-21T09:30:06.861404Z","iopub.status.idle":"2023-07-21T09:30:20.437138Z","shell.execute_reply.started":"2023-07-21T09:30:06.861365Z","shell.execute_reply":"2023-07-21T09:30:20.436088Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(df['term'].value_counts() )\nprint( df['term'].value_counts().describe() )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:20.43884Z","iopub.execute_input":"2023-07-21T09:30:20.439262Z","iopub.status.idle":"2023-07-21T09:30:22.526214Z","shell.execute_reply.started":"2023-07-21T09:30:20.439217Z","shell.execute_reply":"2023-07-21T09:30:22.524768Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print( (df['term'].value_counts() < 200 ).sum() )\nprint( (df['term'].value_counts() < 500 ).sum() )\nprint( (df['term'].value_counts() < 1000 ).sum() )\nprint( (df['term'].value_counts() < 5000 ).sum() )\n\nplt.figure(figsize = (20,4))\nplt.plot(np.log10(1 + v.head(500)),'*-')\nplt.title('Log10 count terms with give number of labels - show terms with least possible number of labels - 1,2,3...', fontsize = 20 )\nplt.xlabel('number of labels', fontsize = 20 )\nplt.grid()\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:22.527712Z","iopub.execute_input":"2023-07-21T09:30:22.528086Z","iopub.status.idle":"2023-07-21T09:30:26.870408Z","shell.execute_reply.started":"2023-07-21T09:30:22.528047Z","shell.execute_reply":"2023-07-21T09:30:26.868954Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"!pip install obonet","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:26.872238Z","iopub.execute_input":"2023-07-21T09:30:26.872602Z","iopub.status.idle":"2023-07-21T09:30:41.260027Z","shell.execute_reply.started":"2023-07-21T09:30:26.872568Z","shell.execute_reply":"2023-07-21T09:30:41.258511Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import obonet\nclass CFG:\n    train_go_obo_path: str = \"/kaggle/input/cafa-5-protein-function-prediction/Train/go-basic.obo\"\n    train_seq_fasta_path: str = \"/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta\"\n    train_terms_path: str = \"/kaggle/input/cafa-5-protein-function-prediction/Train/train_terms.tsv\"\n    train_taxonomy_path: str = \"/kaggle/input/cafa-5-protein-function-prediction/Train/train_taxonomy.tsv\"\n    train_ia_path: str = \"/kaggle/input/cafa-5-protein-function-prediction/IA.txt\"\ngraph = obonet.read_obo(CFG.train_go_obo_path)\nprint(f\"Number of nodes: {len(graph)}\")\nprint(f\"Number of edges: {graph.number_of_edges()}\")","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:30:41.262271Z","iopub.execute_input":"2023-07-21T09:30:41.262757Z","iopub.status.idle":"2023-07-21T09:31:01.781475Z","shell.execute_reply.started":"2023-07-21T09:30:41.262704Z","shell.execute_reply":"2023-07-21T09:31:01.780039Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"term = \"GO:0034655\"\ngraph.nodes[term]","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:01.783275Z","iopub.execute_input":"2023-07-21T09:31:01.783698Z","iopub.status.idle":"2023-07-21T09:31:01.792919Z","shell.execute_reply.started":"2023-07-21T09:31:01.783656Z","shell.execute_reply":"2023-07-21T09:31:01.791292Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Top and least frequent labels ","metadata":{}},{"cell_type":"code","source":"# subontology_roots = {'BPO':'GO:0008150', #  Top2 frequent \n#                      'CCO':'GO:0005575', # Top frequent \n#                      'MFO':'GO:0003674'} # Top4 frequent\n\nprint('Top frequent labels:')\nv = df['term'].value_counts()\nl = v.head(30).index \nfor term in l:\n#     print(term)\n    print(term, v[term], graph.nodes[term]['name'] )\n    \nprint()\nprint('Least frequent labels:')\nl = df['term'].value_counts().tail(10).index \nfor term in l:\n#     print(term)\n    print(term, graph.nodes[term]['name'] )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:01.796859Z","iopub.execute_input":"2023-07-21T09:31:01.797372Z","iopub.status.idle":"2023-07-21T09:31:03.853693Z","shell.execute_reply.started":"2023-07-21T09:31:01.797328Z","shell.execute_reply":"2023-07-21T09:31:03.852329Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Targets and weights ","metadata":{}},{"cell_type":"code","source":"pd.set_option('display.max_colwidth', 200)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:03.856087Z","iopub.execute_input":"2023-07-21T09:31:03.856504Z","iopub.status.idle":"2023-07-21T09:31:03.861902Z","shell.execute_reply.started":"2023-07-21T09:31:03.856462Z","shell.execute_reply":"2023-07-21T09:31:03.860494Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfn = '/kaggle/input/cafa-5-protein-function-prediction/IA.txt'\ndf_ia = pd.read_csv(fn , sep = '\\t', header = None, index_col = 0,  )\ndf_ia.columns = ['Weight']\nprint(df_ia.shape)\ndisplay( df_ia.head(3) )\n\n\ndd  = df['term'].value_counts().to_frame()\ndd.columns = ['count']\ndd = dd.join(df_ia, how = 'outer').sort_values('count',ascending = False).fillna(0)\ndd['count*Weight'] = dd['count']*dd['Weight'].round(2)\nl = [graph.nodes[term]['name'] for term in dd.index ]\ndd['name'] = l\n_t = {'cellular_component':'CCO', 'biological_process':'BPO','molecular_function':'MFO'}\nl = [_t[graph.nodes[term]['namespace']] for term in dd.index ]\ndd['namespace'] = l\nl = [graph.nodes[term]['def'] for term in dd.index ]\ndd['def'] = l\n\ndisplay(dd.head(50))\n\nN = len(dd)# 20_000\nv2 = np.log10(1+dd['count']).values\nl = [dd['Weight'].head(N).max() for N in range(dd.shape[0])]\nplt.figure(figsize = (20,3) )\nplt.plot(l[:N], label = 'cum weights')\nplt.plot(v2[:N], label = 'count' )\nplt.legend()\nplt.title('Cummulative max of weight ordered by frequency')\nplt.grid()\nplt.show()\nl = [dd['Weight'].head(N).mean() for N in range(dd.shape[0])]\nplt.figure(figsize = (20,3) )\nplt.plot(l[:N], label = 'cum weights')\nplt.plot(v2[:N], label = 'count' )\nplt.legend()\nplt.title('Cummulative mean of weight ordered by frequency')\nplt.grid()\nplt.show()\nl = [dd['Weight'].head(N).median() for N in range(dd.shape[0])]\nplt.figure(figsize = (20,3) )\nplt.plot(l[:N], label = 'cum weights')\nplt.plot(v2[:N], label = 'count' )\nplt.legend()\nplt.title('Cummulative median of weight ordered by frequency')\nplt.grid()\nplt.show()\n\nprint(dd.shape)\ndd.to_csv('df_terms_counts_weights_names_etc.csv')\ndisplay(dd)\n\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:03.86378Z","iopub.execute_input":"2023-07-21T09:31:03.864265Z","iopub.status.idle":"2023-07-21T09:31:50.073051Z","shell.execute_reply.started":"2023-07-21T09:31:03.86421Z","shell.execute_reply":"2023-07-21T09:31:50.071437Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Order by count * Weight","metadata":{}},{"cell_type":"code","source":"dd['count*Weight'] = dd['count']*dd['Weight'].round(2)\nplt.figure(figsize = (20,4))\nplt.plot(dd.sort_values('count*Weight', ascending = False)['count*Weight'].head(500).values)\nplt.grid()\nplt.show()\nplt.figure(figsize = (20,4))\nplt.plot(dd.sort_values('count*Weight', ascending = False)['count*Weight'].head(40_000).cumsum().values)\nplt.grid()\nplt.show()\n\ndisplay(dd.sort_values('count*Weight', ascending = False).head(50) )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:50.07478Z","iopub.execute_input":"2023-07-21T09:31:50.075887Z","iopub.status.idle":"2023-07-21T09:31:50.691805Z","shell.execute_reply.started":"2023-07-21T09:31:50.075845Z","shell.execute_reply":"2023-07-21T09:31:50.690348Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"N= 30 \nfor t in [15_000,10_000,5_000, 1000,100]:\n    m = dd['count'] > t\n    print(); print(); print(t,m.sum() )\n    display( dd[m].sort_values('Weight', ascending = False ).head(N) )\n\n\ndisplay( dd.head(1500).sort_values('Weight', ascending = False ).head(N) )\n    ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:50.693359Z","iopub.execute_input":"2023-07-21T09:31:50.693886Z","iopub.status.idle":"2023-07-21T09:31:50.861002Z","shell.execute_reply.started":"2023-07-21T09:31:50.693832Z","shell.execute_reply":"2023-07-21T09:31:50.859323Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print(df['aspect'].value_counts())\nprint()\nv = df['term'].value_counts()\n\nfor t in [1,10]:\n    m = v <= t\n    print()\n    print(t,m.sum())\n    l = []\n    for term in v[m].index:\n        l.append( graph.nodes[term]['namespace']  )\n    print(pd.Series(l).value_counts() )\n\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:50.862886Z","iopub.execute_input":"2023-07-21T09:31:50.863402Z","iopub.status.idle":"2023-07-21T09:31:52.927229Z","shell.execute_reply.started":"2023-07-21T09:31:50.863349Z","shell.execute_reply":"2023-07-21T09:31:52.925729Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Cell cycle terms","metadata":{}},{"cell_type":"code","source":"m = np.array([ 'cell cycle' in t for t in dd['name']])\nprint(m.sum())\ndisplay(dd[m].iloc[:50])\ndisplay(dd[m].iloc[50:100])\ndisplay(dd[m].iloc[100:])","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:52.929224Z","iopub.execute_input":"2023-07-21T09:31:52.929675Z","iopub.status.idle":"2023-07-21T09:31:53.063831Z","shell.execute_reply.started":"2023-07-21T09:31:52.929614Z","shell.execute_reply":"2023-07-21T09:31:53.062482Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Sequences not starting from \"M\" letter (Methionine) \n\n\nProtein synthesis is believed to be initiated with the amino acid methionine because the AUG translation initiation codon of mRNAs is recognized by the anticodon of initiator methionine transfer RNA. A group of positive-stranded RNA viruses of insects, however, lacks an AUG translation initiation codon for their capsid protein gene, which is located at the downstream part of the genome. The capsid protein of one of these viruses, Plautia stali intestine virus, is synthesized by internal ribosome entry site-mediated translation. Here we report that methionine is not the initiating amino acid in the translation of the capsid protein in this virus. Its translation is initiated with glutamine encoded by a CAA codon that is the first codon of the capsid-coding region. \n\nhttps://doi.org/10.1073/pnas.010426997\n","metadata":{}},{"cell_type":"markdown","source":"## In train ( look: M is not the first)","metadata":{}},{"cell_type":"code","source":"%%time\nprint('Train data:')\nfile_fasta = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'#  '/kaggle/input/biopython-genbank/NC_005816.fna'\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nc = 0\nl = []; l2 = []; l3_train = []; l4_train= [];\nfor i,seq in list(enumerate(sequences)):\n    if ( 'M' != str(seq.seq)[0]):\n        l.append(seq.description); l2.append(len(seq.seq)); l3_train.append(str(seq.seq)); l4_train.append(seq.id);\n    if (c<10) and ( 'M' != str(seq.seq)[0]):\n        c += 1\n        print(seq.id, len(seq.seq), '    ', seq.seq[:30]);         print(seq.description)\nprint(); print('count:', len(l)); print('Mean len: %.1f, median %.1f'%(np.mean(l2), np.median(l2)) ); l_train = l.copy()\n\nl2 = []\nfor s in l:\n    if ('peptide' not in s.lower() ) and ( 'fragment' not in s.lower() ):\n        l2.append(s)\nprint(); print('count:', len(l2))   \nl2 = []\nfor s in l:\n    if ('virus' not in s.lower() ) :\n        l2.append(s)\nprint(); print('count:', len(l2))        ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:53.065867Z","iopub.execute_input":"2023-07-21T09:31:53.06668Z","iopub.status.idle":"2023-07-21T09:31:58.045446Z","shell.execute_reply.started":"2023-07-21T09:31:53.066599Z","shell.execute_reply":"2023-07-21T09:31:58.04396Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## In test  ( look: M is not the first)","metadata":{}},{"cell_type":"code","source":"%%time\nprint('Test data: \\n')\nfile_fasta2 = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences2 = SeqIO.parse(file_fasta2, \"fasta\")\nc = 0\nl = []; l2 = []; l3_test = []; l4_test = [] \nfor i,seq in list(enumerate(sequences2)):\n    if ( 'M' != str(seq.seq)[0]):\n        l.append(seq.description); l2.append(len(seq)); l3_test.append(str(seq.seq));  l4_test.append(seq.id);\n    if (c<10) and ( 'M' != str(seq.seq)[0]):\n        c += 1\n        print(seq.id, len(seq.seq), '    ', seq.seq[:30]);         print(seq.description)\nprint(); print('count:', len(l)); print('Mean len: %.1f, median %.1f'%(np.mean(l2), np.median(l2)) ); l_test = l.copy()\nl2 = []\nfor s in l:\n    if ('peptide' not in s.lower() ) and ( 'fragment' not in s.lower() ):\n        l2.append(s)\nprint(); print('count:', len(l2))   \nl2 = []\nfor s in l:\n    if ('virus' not in s.lower() ) :\n        l2.append(s)\nprint(); print('count:', len(l2))        ","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:31:58.047772Z","iopub.execute_input":"2023-07-21T09:31:58.048293Z","iopub.status.idle":"2023-07-21T09:32:02.714065Z","shell.execute_reply.started":"2023-07-21T09:31:58.04824Z","shell.execute_reply":"2023-07-21T09:32:02.712697Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Intersection train test","metadata":{}},{"cell_type":"code","source":"s = set(l3_train) & set(l3_test)\nprint(len(s))\ns = set(l4_train) & set(l4_test)\nprint(len(s))\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:02.715619Z","iopub.execute_input":"2023-07-21T09:32:02.716021Z","iopub.status.idle":"2023-07-21T09:32:02.729322Z","shell.execute_reply.started":"2023-07-21T09:32:02.715985Z","shell.execute_reply":"2023-07-21T09:32:02.727729Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# Duplicated sequences ","metadata":{}},{"cell_type":"code","source":"%%time \n# Train\nfile_fasta = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'#  '/kagg\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nlist_ids = [seq.id for seq in sequences]\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nlist_lens = [len(seq.seq) for seq in sequences]\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nlist_seqs = [str(seq.seq) for seq in sequences]\nprint('Mean len: %.1f, median %.1f'%(np.mean(list_lens), np.median(list_lens)) )\n\n# Test\nfile_fasta2 = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences = SeqIO.parse(file_fasta2, \"fasta\")\nlist_ids2 = [seq.id for seq in sequences]\nsequences = SeqIO.parse(file_fasta2, \"fasta\")\nlist_lens2 = [len(seq.seq) for seq in sequences]\nsequences = SeqIO.parse(file_fasta2, \"fasta\")\nlist_seqs2 = [str(seq.seq) for seq in sequences]\nprint('Mean len: %.1f, median %.1f'%(np.mean(list_lens2), np.median(list_lens2)) )\n\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:02.731059Z","iopub.execute_input":"2023-07-21T09:32:02.731487Z","iopub.status.idle":"2023-07-21T09:32:18.589982Z","shell.execute_reply.started":"2023-07-21T09:32:02.731447Z","shell.execute_reply":"2023-07-21T09:32:18.588281Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Duplicated sequences in train","metadata":{}},{"cell_type":"code","source":"%%time\nv = pd.Series(list_seqs).value_counts()# .head(10)\nprint('total duplicated sequences:', v[v!=1].sum())\nplt.figure(figsize = (20,6) )\nplt.plot(v.values[:100],'*-')\nplt.title('Duplicate sequences statistics Train', fontsize = 20 )\nplt.ylabel('number of duplicates')\nplt.xlabel('index')\nplt.grid()\nplt.show()\nfor t in [2,3,4,5,6]:\n    print('number of proteins with ',t, 'dubplicates:', (v==t).sum() )\nprint((v>1).sum(), list(v[v>1].values[:350]) )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:18.593785Z","iopub.execute_input":"2023-07-21T09:32:18.59469Z","iopub.status.idle":"2023-07-21T09:32:19.143968Z","shell.execute_reply.started":"2023-07-21T09:32:18.59462Z","shell.execute_reply":"2023-07-21T09:32:19.142739Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### examples of duplicates ","metadata":{}},{"cell_type":"code","source":"%%time\nprint(v.iloc[100:101])\nseq1 = 'MENRNTFSWVKEQMTRSISVSIMIYVITRTSISNAYPIFAQQGYENPREATGRIVCANCHLANKPVDIEVPQAVLPDTVFEAVLRIPYDMQLKQVLANGKKGGLNVGAVLILPEGFELAPPDRISPELKEKIGNLSFQSYRPNKKNILVIGPVPGKKYSEIVFPILSPDPAMKKDVHFLKYPIYVGGNRGRGQIYPDGSKSNNTVYNATSTGVVRKILRKEKGGYEISIVDASDGRQVIDLIPPGPELLVSEGESIKLDQPLTSNPNVGGFGQGDAEIVLQDPLRVQGLLFFFASVILAQVFLVLKKKQFEKVQLYEMNF'\n\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nfor i,seq in list(enumerate(sequences)):\n    if seq1 == seq.seq:\n        print(i, seq.id, len(seq.seq))\n        print(seq.seq[:50])\n        print(seq.description)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:19.145917Z","iopub.execute_input":"2023-07-21T09:32:19.14677Z","iopub.status.idle":"2023-07-21T09:32:24.282932Z","shell.execute_reply.started":"2023-07-21T09:32:19.146728Z","shell.execute_reply":"2023-07-21T09:32:24.281435Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### top duplicated case ","metadata":{}},{"cell_type":"code","source":"%%time\nprint( v.iloc[:1])\nseq1 = 'MTQSNPNEQNVELNRTSLYWGLLLIFVLAVLFSNYFFN'\n\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nfor i,seq in list(enumerate(sequences)):\n    if seq1 == seq.seq:\n        print(i, seq.id, len(seq.seq))\n        print(seq.seq[:50])\n        print(seq.description)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:24.285066Z","iopub.execute_input":"2023-07-21T09:32:24.286452Z","iopub.status.idle":"2023-07-21T09:32:28.833549Z","shell.execute_reply.started":"2023-07-21T09:32:24.286395Z","shell.execute_reply":"2023-07-21T09:32:28.832098Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Duplicated sequences in test","metadata":{}},{"cell_type":"code","source":"%%time\nv = pd.Series(list_seqs2).value_counts()# .head(10)\nprint('total duplicated sequences:', v[v!=1].sum())\nplt.figure(figsize = (20,6) )\nplt.plot(v.values[:100],'*-')\nplt.title('Duplicate sequences statistics Test', fontsize = 20 )\nplt.ylabel('number of duplicates')\nplt.xlabel('index')\nplt.grid()\nplt.show()\nfor t in [2,3,4,5,6]:\n    print('number of proteins with ',t, 'dubplicates:', (v==t).sum() )\nprint((v>1).sum(), list(v[v>1].values[:350]) )\n\nprint()\nprint('Info on top duplicated protein:')\nprint( v.iloc[:1])\nseq1 = 'MSHQLTFADSEFSSKRRQTRKEIFLSRMEQILPWQNMVEVIEPFYPKAGNGRRPYPLETMLRIHCMQHWYNLSDGAMEDALYEIASMRLFARLSLDSALPDRTTIMNFRHLLEQHQLARQLFKTINRWLAEAGVMMTQGTLVDATIIEAPSSTKNKEQQRDPEMHQTKKGNQWHFGMKAHIGVDAKSGLTHSLVTTAANEHDLNQLGNLLHGEEQFVSADAGYQGAPQREELAEVDVDWLIAERPGKVRTLKQHPRKNKTAINIEYMKASIRARVEHPFRIIKRQFGFVKARYKGLLKNDNQLAMLFTLANLFRADQMIRQWERSH'\nsequences = SeqIO.parse(file_fasta2, \"fasta\")\nfor i,seq in list(enumerate(sequences)):\n    if seq1 == seq.seq:\n        print(i, seq.id, len(seq.seq), 'In train:', str(seq.seq) in list_seqs)\n#         print(seq.seq[:50])\n#         print(seq.description)\nprint()\nprint('Info on top-2 duplicated protein:')\nprint( v.iloc[1:2])\nseq1 = v.index[1]\nsequences = SeqIO.parse(file_fasta2, \"fasta\")\nfor i,seq in list(enumerate(sequences)):\n    if seq1 == seq.seq:\n        print(i, seq.id, len(seq.seq), 'In train:', str(seq.seq) in list_seqs)\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:28.837914Z","iopub.execute_input":"2023-07-21T09:32:28.838319Z","iopub.status.idle":"2023-07-21T09:32:38.98839Z","shell.execute_reply.started":"2023-07-21T09:32:28.838282Z","shell.execute_reply":"2023-07-21T09:32:38.987017Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"# ","metadata":{}},{"cell_type":"code","source":"%%time\nfn = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_taxonomy.tsv'\ndf = pd.read_csv(fn, sep = '\\t')\ndisplay(df.head(2))\nv = df['taxonomyID'].value_counts()\nv","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:38.990394Z","iopub.execute_input":"2023-07-21T09:32:38.99094Z","iopub.status.idle":"2023-07-21T09:32:39.101968Z","shell.execute_reply.started":"2023-07-21T09:32:38.990881Z","shell.execute_reply":"2023-07-21T09:32:39.100466Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from Bio import Entrez\n\n# Provide your email address to the Entrez module (required by NCBI)\nEntrez.email = \"your_email@example.com\"\n\ndef get_taxonomic_name(taxon_id):\n    handle = Entrez.efetch(db=\"taxonomy\", id=str(taxon_id), retmode=\"xml\")\n    records = Entrez.read(handle)\n    handle.close()\n\n    for record in records:\n        name = record.get(\"ScientificName\", None)\n        if name:\n            return name\n\n    return None\n\n# Example usage\ntaxon_id = 9606  # Taxon ID for Homo sapiens (human)\ntaxonomic_name = get_taxonomic_name(taxon_id)\n# print(taxonomic_name)\n\ndf2 = pd.DataFrame()\ndf2['taxonomyID'] = v.index\ndf2['count'] = v.values\ndf2['name'] = ''\nfor i,k in enumerate(v.index):\n    taxonomic_name = get_taxonomic_name(k)\n    df2.loc[i,'name'] = taxonomic_name\n    if i%5 == 0: print(i,taxonomic_name)","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:32:39.103489Z","iopub.execute_input":"2023-07-21T09:32:39.104534Z","iopub.status.idle":"2023-07-21T09:33:24.357254Z","shell.execute_reply.started":"2023-07-21T09:32:39.104489Z","shell.execute_reply":"2023-07-21T09:33:24.354235Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"df2","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:33:24.358854Z","iopub.status.idle":"2023-07-21T09:33:24.360221Z","shell.execute_reply.started":"2023-07-21T09:33:24.359923Z","shell.execute_reply":"2023-07-21T09:33:24.359959Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%time\nfile_fasta = '/kaggle/input/cafa-5-protein-function-prediction/Train/train_sequences.fasta'#  '/kaggle/input/biopython-genbank/NC_005816.fna'\nsequences = SeqIO.parse(file_fasta, \"fasta\")\nlist_ids = [seq.id for seq in sequences]\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:36:16.112601Z","iopub.execute_input":"2023-07-21T09:36:16.113123Z","iopub.status.idle":"2023-07-21T09:36:18.737029Z","shell.execute_reply.started":"2023-07-21T09:36:16.113083Z","shell.execute_reply":"2023-07-21T09:36:18.735558Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"file_fasta2 = '/kaggle/input/cafa-5-protein-function-prediction/Test (Targets)/testsuperset.fasta'\nsequences = SeqIO.parse(file_fasta2, \"fasta\")\nlist_ids2 = [seq.id for seq in sequences]\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:36:18.739232Z","iopub.execute_input":"2023-07-21T09:36:18.739726Z","iopub.status.idle":"2023-07-21T09:36:21.308571Z","shell.execute_reply.started":"2023-07-21T09:36:18.739687Z","shell.execute_reply":"2023-07-21T09:36:21.307273Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"for i in  ['P04637', 'O15350', 'P06493', 'Q9Y3S1', 'Q96J92', 'Q9BYP7' , 'Q9H4A3', 'Q99592']: # TP53, TP73 , CDK1, WNK2 WNK4 WNK3 WNK1 ZBTB18\n    print(i, i in list_ids, i in list_ids2 )\n","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:36:21.310028Z","iopub.execute_input":"2023-07-21T09:36:21.310391Z","iopub.status.idle":"2023-07-21T09:36:21.384269Z","shell.execute_reply.started":"2023-07-21T09:36:21.310355Z","shell.execute_reply":"2023-07-21T09:36:21.382837Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"print('%.1f seconds passed total '%(time.time()-t0start) )","metadata":{"execution":{"iopub.status.busy":"2023-07-21T09:36:31.023409Z","iopub.execute_input":"2023-07-21T09:36:31.024779Z","iopub.status.idle":"2023-07-21T09:36:31.030597Z","shell.execute_reply.started":"2023-07-21T09:36:31.024727Z","shell.execute_reply":"2023-07-21T09:36:31.029527Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}