{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"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":"gpu","dataSources":[{"sourceId":116062,"databundleVersionId":14084779,"sourceType":"competition"}],"dockerImageVersionId":31193,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":true}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"#Imports\n\nimport os, glob, gzip, re\nimport pandas as pd, numpy as np\nfrom sklearn.feature_extraction.text import CountVectorizer\nfrom sklearn.linear_model import LogisticRegression\nfrom sklearn.multiclass import OneVsRestClassifier\nfrom sklearn.preprocessing import MultiLabelBinarizer\nfrom sklearn.model_selection import train_test_split\nfrom sklearn.metrics import f1_score\nfrom sklearn.feature_extraction.text import HashingVectorizer\nfrom sklearn.linear_model import SGDClassifier\nfrom scipy import sparse\nimport math\n\n!pip install -q biopython\n\n# test import immediately\ntry:\n    from Bio import SeqIO\n    print('Biopython installed and import OK:', SeqIO.__name__)\nexcept Exception as e:\n    print('Import failed — you may need to restart the kernel. Error:', e) ","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:41.518425Z","iopub.execute_input":"2025-11-16T18:50:41.519059Z","iopub.status.idle":"2025-11-16T18:50:47.916711Z","shell.execute_reply.started":"2025-11-16T18:50:41.519033Z","shell.execute_reply":"2025-11-16T18:50:47.915783Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#File Detection\n\nINPUT_ROOT = '/kaggle/input'\nall_files = glob.glob(os.path.join(INPUT_ROOT, '**', '*'), recursive=True)\n\nfasta_train = None\nfasta_test = None\nterms_tsv = None\n\nfor p in all_files:\n    low = p.lower()\n    if low.endswith('.fasta') or low.endswith('.fa') or low.endswith('.fasta.gz') or low.endswith('.fa.gz'):\n        if 'train' in low and fasta_train is None:\n            fasta_train = p\n        elif 'test' in low and fasta_test is None:\n            fasta_test = p\n        elif fasta_train is None:\n            fasta_train = p\n    if (low.endswith('.tsv') or low.endswith('.csv')) and terms_tsv is None:\n        if 'term' in low or 'label' in low or 'annotation' in low or 'train_terms' in low:\n            terms_tsv = p\n\ncomp_folder = os.path.join(INPUT_ROOT, 'cafa-6-protein-function-prediction')\nif os.path.exists(comp_folder):\n    if fasta_train is None and os.path.exists(os.path.join(comp_folder, 'train_sequences.fasta')):\n        fasta_train = os.path.join(comp_folder, 'train_sequences.fasta')\n    if fasta_test is None and os.path.exists(os.path.join(comp_folder, 'test_sequences.fasta')):\n        fasta_test = os.path.join(comp_folder, 'test_sequences.fasta')\n    if terms_tsv is None and os.path.exists(os.path.join(comp_folder, 'train_terms.tsv')):\n        terms_tsv = os.path.join(comp_folder, 'train_terms.tsv')\n\nif not fasta_train or not fasta_test or not terms_tsv:\n    print('Auto-detection failed. Found:')\n    print(' fasta_train =', fasta_train)\n    print(' fasta_test  =', fasta_test)\n    print(' terms_tsv   =', terms_tsv)\n    raise SystemExit('Attach the CAFA-6 dataset: Add data -> cafa-6-protein-function-prediction')\n\nprint('Using files:')\nprint(' train fasta:', fasta_train)\nprint(' test  fasta:', fasta_test)\nprint(' terms tsv :', terms_tsv)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:47.918413Z","iopub.execute_input":"2025-11-16T18:50:47.919027Z","iopub.status.idle":"2025-11-16T18:50:47.937269Z","shell.execute_reply.started":"2025-11-16T18:50:47.918998Z","shell.execute_reply":"2025-11-16T18:50:47.936662Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#Defining Helpers\n\ndef read_fasta_ids(path, limit=None):\n    open_fn = gzip.open if str(path).endswith('.gz') else open\n    ids = []\n    with open_fn(path, 'rt') as handle:\n        for rec in SeqIO.parse(handle, 'fasta'):\n            ids.append(rec.id)\n            if limit and len(ids) >= limit:\n                break\n    return ids\n\ndef read_fasta_dict(path):\n    open_fn = gzip.open if str(path).endswith('.gz') else open\n    seqs = {}\n    with open_fn(path, 'rt') as handle:\n        for rec in SeqIO.parse(handle, 'fasta'):\n            seqs[rec.id] = str(rec.seq)\n    return seqs","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:47.937947Z","iopub.execute_input":"2025-11-16T18:50:47.938191Z","iopub.status.idle":"2025-11-16T18:50:47.943619Z","shell.execute_reply.started":"2025-11-16T18:50:47.938169Z","shell.execute_reply":"2025-11-16T18:50:47.942992Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#previewing and debudding\n\ntrain_sample_ids = read_fasta_ids(fasta_train, limit=20)\ntest_sample_ids = read_fasta_ids(fasta_test, limit=20)\nprint('\\nExample train FASTA ids (20):', train_sample_ids[:20])\nprint('Example test  FASTA ids (20):', test_sample_ids[:20])","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:47.945014Z","iopub.execute_input":"2025-11-16T18:50:47.945257Z","iopub.status.idle":"2025-11-16T18:50:47.970619Z","shell.execute_reply.started":"2025-11-16T18:50:47.945241Z","shell.execute_reply":"2025-11-16T18:50:47.969951Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#read and detect terms_tsv format\n\nopen_fn = gzip.open if str(terms_tsv).endswith('.gz') else open\nprint('\\n--- first 12 raw lines of terms_tsv ---')\nwith open_fn(terms_tsv, 'rt') as f:\n    for i, line in enumerate(f):\n        print(i+1, line.strip())\n        if i >= 11:\n            break","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:47.971259Z","iopub.execute_input":"2025-11-16T18:50:47.971456Z","iopub.status.idle":"2025-11-16T18:50:47.978385Z","shell.execute_reply.started":"2025-11-16T18:50:47.971440Z","shell.execute_reply":"2025-11-16T18:50:47.977761Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#trying with more rows \n\ndf = pd.read_csv(terms_tsv, sep='\\t', header=None, dtype=str, engine='python')\nprint('\\nRead shape (terms file):', df.shape)\nprint('First rows:')\nprint(df.head(6))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:47.979226Z","iopub.execute_input":"2025-11-16T18:50:47.979850Z","iopub.status.idle":"2025-11-16T18:50:49.548788Z","shell.execute_reply.started":"2025-11-16T18:50:47.979830Z","shell.execute_reply":"2025-11-16T18:50:49.547988Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"# heuristic- check for GO or similar\n\ncol0_is_go = df[0].astype(str).str.match(r'^GO:\\\\d{7}').sum() if 0 in df.columns else 0\ncol1_is_go = df[1].astype(str).str.match(r'^GO:\\\\d{7}').sum() if 1 in df.columns else 0\nprint(f'col0 GO-like count = {col0_is_go}, col1 GO-like count = {col1_is_go}')\n\nif col0_is_go > col1_is_go:\n    print('Detected first column contains GO terms. Interpreting file as (term_id, protein_id) -> swapping.')\n    terms_df = df.rename(columns={0:'term_id', 1:'protein_id'})[['protein_id','term_id']]\nelse:\n    print('Detected first column contains protein ids. Using (protein_id, term_id).')\n    terms_df = df.rename(columns={0:'protein_id', 1:'term_id'})[['protein_id','term_id']]\n\nprint('\\nterms_df sample:')\nprint(terms_df.head(6))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:49.549685Z","iopub.execute_input":"2025-11-16T18:50:49.549971Z","iopub.status.idle":"2025-11-16T18:50:50.067413Z","shell.execute_reply.started":"2025-11-16T18:50:49.549947Z","shell.execute_reply":"2025-11-16T18:50:50.066559Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#Normalization \n\ndef norm_variants(s):\n    s = '' if s is None else str(s)\n    out = []\n    out.append(s)\n    if '|' in s:\n        parts = s.split('|')\n        for p in parts:\n            if p: out.append(p)\n        out.append(parts[-1])\n        if len(parts) > 1: out.append(parts[1])\n    out.append(s.split()[0])\n    if '.' in s: out.append(s.split('.')[0])\n    out.append(re.sub('[^A-Za-z0-9_\\\\-]', '', s))\n    for t in re.split('[\\\\|\\\\s]', s):\n        t = t.strip()\n        if 4 <= len(t) <= 12:\n            out.append(t)\n    uniq = []\n    for v in out:\n        if v and v not in uniq: uniq.append(v)\n    return uniq\n\nlabel_ids = pd.unique(terms_df['protein_id'])\nnorm2labels = {}\nfor lid in label_ids:\n    for nv in norm_variants(lid):\n        norm2labels.setdefault(nv, set()).add(lid)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:50.068318Z","iopub.execute_input":"2025-11-16T18:50:50.068615Z","iopub.status.idle":"2025-11-16T18:50:50.453306Z","shell.execute_reply.started":"2025-11-16T18:50:50.068589Z","shell.execute_reply":"2025-11-16T18:50:50.452506Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"print('\\nSample match count (train sample):', sum(1 for tid in train_sample_ids if any(c in norm2labels for c in norm_variants(tid))), '/', len(train_sample_ids))\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:50.454420Z","iopub.execute_input":"2025-11-16T18:50:50.454994Z","iopub.status.idle":"2025-11-16T18:50:50.459444Z","shell.execute_reply.started":"2025-11-16T18:50:50.454970Z","shell.execute_reply":"2025-11-16T18:50:50.458729Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#Map train (fasta ids -> lable ids)\n\ntrain_all_ids = read_fasta_ids(fasta_train, limit=None)\nmapped_train_to_label = {}\nunmatched = []\nfor tid in train_all_ids:\n    found = False\n    for cand in norm_variants(tid):\n        if cand in norm2labels:\n            mapped_train_to_label[tid] = list(norm2labels[cand])[0]\n            found = True\n            break\n    if not found:\n        unmatched.append(tid)\n\nprint('Total train sequences:', len(train_all_ids))\nprint('Mapped train->label count:', len(mapped_train_to_label))\nprint('Unmatched train sequences:', len(unmatched))\nprint('Example unmatched (10):', unmatched[:10])\n\nif len(mapped_train_to_label) == 0:\n    raise SystemExit('No mapping between fasta ids and label ids found automatically. Paste the first raw lines above here for further rule craft.')\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:50.461320Z","iopub.execute_input":"2025-11-16T18:50:50.461683Z","iopub.status.idle":"2025-11-16T18:50:51.634031Z","shell.execute_reply.started":"2025-11-16T18:50:50.461666Z","shell.execute_reply":"2025-11-16T18:50:51.633299Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#mapping terms_df\n\nlabel_to_fasta = {}\nfor f_id, l_id in mapped_train_to_label.items():\n    label_to_fasta.setdefault(l_id, set()).add(f_id)\n\ndef map_label_row(lid):\n    if lid in label_to_fasta:\n        return list(label_to_fasta[lid])[0]\n    for nv in norm_variants(lid):\n        if nv in label_to_fasta:\n            return list(label_to_fasta[nv])[0]\n    return None\n\nterms_df['mapped_fasta_id'] = terms_df['protein_id'].apply(map_label_row)\nmapped_count = int(terms_df['mapped_fasta_id'].notnull().sum())\nprint('Mapped label rows to fasta ids:', mapped_count, '/', len(terms_df))\n\nmapped_df = terms_df[terms_df['mapped_fasta_id'].notnull()].copy()\ngrouped = mapped_df.groupby('mapped_fasta_id')['term_id'].apply(list).reset_index().rename(columns={'mapped_fasta_id':'protein_id'})\nprint('Grouped unique proteins with labels (after mapping):', len(grouped))","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:51.634708Z","iopub.execute_input":"2025-11-16T18:50:51.634885Z","iopub.status.idle":"2025-11-16T18:50:53.480419Z","shell.execute_reply.started":"2025-11-16T18:50:51.634871Z","shell.execute_reply":"2025-11-16T18:50:53.479658Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#loading seq dict\n\ntrain_seqs = read_fasta_dict(fasta_train)\ntest_seqs = read_fasta_dict(fasta_test)\n\ngrouped = grouped[grouped['protein_id'].isin(train_seqs.keys())].reset_index(drop=True)\nprint('Grouped after filtering to available sequences:', len(grouped))\nif len(grouped) == 0:\n    raise SystemExit('No labelled proteins remain after mapping and filtering. Inspect unmatched examples above.')\n\ngrouped['sequence'] = grouped['protein_id'].map(train_seqs)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:53.481247Z","iopub.execute_input":"2025-11-16T18:50:53.481499Z","iopub.status.idle":"2025-11-16T18:50:55.417287Z","shell.execute_reply.started":"2025-11-16T18:50:53.481474Z","shell.execute_reply":"2025-11-16T18:50:55.416520Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#k-mers\n\nK = 3\ndef kmers(seq, k=K):\n    if not isinstance(seq, str) or len(seq) == 0: return ['']\n    if len(seq) < k: return [seq]\n    return [seq[i:i+k] for i in range(len(seq)-k+1)]\n\ntrain_texts = [' '.join(kmers(s)) for s in grouped['sequence']]\n\nvectorizer = HashingVectorizer(n_features=2**15, alternate_sign=False, token_pattern=r'[^\\\\s]+')  # 32k features\nX_full = vectorizer.transform(train_texts)\nif 'mlb' not in globals():\n    from sklearn.preprocessing import MultiLabelBinarizer\n    mlb = MultiLabelBinarizer(sparse_output=True)\nY_sparse = mlb.fit_transform(grouped['term_id'])\nprint('X_full shape:', X_full.shape, 'Y shape:', Y_sparse.shape)","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:50:55.418125Z","iopub.execute_input":"2025-11-16T18:50:55.418308Z","iopub.status.idle":"2025-11-16T18:51:05.745300Z","shell.execute_reply.started":"2025-11-16T18:50:55.418294Z","shell.execute_reply":"2025-11-16T18:51:05.744490Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#Training\n\nSAMPLE_FRAC = 0.2             \nn = X_full.shape[0]\nsample_n = max(2000, int(n * SAMPLE_FRAC))\nrng = np.random.RandomState(42)\nsample_idx = rng.choice(n, size=sample_n, replace=False)\nprint('Training on sample:', sample_n, '/', n)\n\nX_sub = X_full[sample_idx]\nY_sub = Y_sparse[sample_idx].toarray()   \n\ncol_sums = Y_sub.sum(axis=0)\nconst_all_zero = np.where(col_sums == 0)[0]\nconst_all_one  = np.where(col_sums == Y_sub.shape[0])[0]\nconst_cols = np.unique(np.concatenate([const_all_zero, const_all_one]))\nprint('Constant columns in sample (will remove):', len(const_cols))\n\nkeep_mask = np.ones(Y_sub.shape[1], dtype=bool)\nkeep_mask[const_cols] = False\nY_sub_reduced = Y_sub[:, keep_mask]\nkept_classes = [c for i,c in enumerate(mlb.classes_) if keep_mask[i]]\nremoved_classes = [c for i,c in enumerate(mlb.classes_) if not keep_mask[i]]\nprint('Kept classes:', len(kept_classes), 'Removed (constant):', len(removed_classes))\n\novr = OneVsRestClassifier(SGDClassifier(loss='log_loss', max_iter=1000, tol=1e-3), n_jobs=1)\novr.fit(X_sub, Y_sub_reduced)\nprint('Trained OVR on sample.')\n\nsaved_state = {\n    'ovr': ovr,\n    'vectorizer': vectorizer,\n    'mlb': mlb,\n    'keep_mask': keep_mask,\n    'const_all_one_idx': const_all_one.tolist()\n}","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2025-11-16T18:51:05.746140Z","iopub.execute_input":"2025-11-16T18:51:05.746381Z","iopub.status.idle":"2025-11-16T19:07:25.573588Z","shell.execute_reply.started":"2025-11-16T18:51:05.746364Z","shell.execute_reply":"2025-11-16T19:07:25.572705Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"#predictions + submissions\n\nimport csv, math, numpy as np\nfrom math import ceil\n\nTHRESH = 0.5\nTOP_K_IF_NONE = 3\nBATCH_SIZE = 2048   \n\nif 'saved_state' not in globals():\n    raise RuntimeError('saved_state not found. Run the training block that sets saved_state before this cell.')\n\ntry:\n    read_fasta_dict\nexcept NameError:\n    import gzip\n    from Bio import SeqIO\n    def read_fasta_dict(path):\n        open_fn = gzip.open if str(path).endswith('.gz') else open\n        seqs = {}\n        with open_fn(path, 'rt') as handle:\n            for rec in SeqIO.parse(handle, 'fasta'):\n                seqs[rec.id] = str(rec.seq)\n        return seqs\n\nif 'test_seqs' not in globals():\n    if 'fasta_test' not in globals():\n        raise RuntimeError('fasta_test path not found. Ensure file detection ran and fasta_test is defined.')\n    print('Loading test sequences (this may take a moment)...')\n    test_seqs = read_fasta_dict(fasta_test)\n    print('Loaded test seqs:', len(test_seqs))\n\novr = saved_state['ovr']\nvectorizer = saved_state['vectorizer']\nmlb = saved_state['mlb']\nkeep_mask = np.array(saved_state['keep_mask'], dtype=bool)\nconst_all_one_idx = set(saved_state.get('const_all_one_idx', []))\n\nn_total = len(mlb.classes_)\nkept_indices = np.where(keep_mask)[0]       \nkept_pos = -np.ones(n_total, dtype=int)\nkept_pos[kept_indices] = np.arange(len(kept_indices))\n\nif 'kmers' not in globals():\n    def kmers(seq, k=3):\n        if not isinstance(seq, str) or len(seq) == 0: return ['']\n        if len(seq) < k: return [seq]\n        return [seq[i:i+k] for i in range(len(seq)-k+1)]\n\ntest_ids = list(test_seqs.keys())\nn_test = len(test_ids)\nn_batches = ceil(n_test / BATCH_SIZE)\nout_path = \"/kaggle/working/submission.tsv\"\n\nwith open(out_path, \"w\", newline=\"\") as f:\n    writer = csv.writer(f, delimiter=\"\\t\")\n    writer.writerow([\"accession\", \"term\"])\n\nprint(f\"Predicting {n_test} proteins in {n_batches} batches (batch_size={BATCH_SIZE})...\")\n\nfor b in range(n_batches):\n    start = b * BATCH_SIZE\n    end = min(n_test, (b + 1) * BATCH_SIZE)\n    batch_ids = test_ids[start:end]\n\n    batch_texts = [' '.join(kmers(test_seqs[pid])) for pid in batch_ids]\n\n    Xb = vectorizer.transform(batch_texts)\n    probs_kept = ovr.predict_proba(Xb)\n\n    if probs_kept.ndim == 1:\n        probs_kept = probs_kept.reshape(-1, 1)\n\n    rows = []\n    for i, pid in enumerate(batch_ids):\n        pk = probs_kept[i]\n\n        local_pass = np.where(pk >= THRESH)[0]\n        if local_pass.size > 0:\n            selected_global = kept_indices[local_pass].tolist()\n        else:\n            top_local = np.argsort(pk)[-TOP_K_IF_NONE:][::-1]\n            selected_global = kept_indices[top_local].tolist()\n\n        for gi in const_all_one_idx:\n            if gi not in selected_global:\n                selected_global.append(gi)\n\n        for gi in selected_global:\n            term_id = mlb.classes_[gi]\n            rows.append((pid, term_id))\n\n    with open(out_path, \"a\", newline=\"\") as f:\n        writer = csv.writer(f, delimiter=\"\\t\")\n        writer.writerows(rows)\n\n    print(f\"Batch {b+1}/{n_batches}  |  Proteins {start}..{end-1}  |  Wrote {len(rows)} rows\")\n\nprint(\"FINAL SUBMISSION SAVED:\", out_path)","metadata":{"trusted":true},"outputs":[],"execution_count":null}]}