{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.12.13","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[],"dockerImageVersionId":28755,"isInternetEnabled":false,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"code","source":"import os,random\nfrom pathlib import Path\nimport numpy as np\nimport pandas as pd\nimport pydicom\nimport matplotlib.pyplot as plt\n\nROOT=Path(\"/kaggle/input/competitions/rsna-knee-abnormality-detection\")\nTRAIN_DIR=ROOT/\"train_series\"\ntrain=pd.read_csv(ROOT/\"train.csv\")\nN_STUDIES=6\nN_SERIES=6\nrandom.seed(1130)\n\nstudies=random.sample(list(train[\"StudyInstanceUID\"].astype(str)),N_STUDIES)\n\ndef read_middle_slice(series_dir):\n    files=sorted(series_dir.glob(\"*.dcm\"))\n    if not files:return None,None\n    ds=pydicom.dcmread(files[len(files)//2],force=True)\n    img=ds.pixel_array.astype(np.float32)\n    slope=float(getattr(ds,\"RescaleSlope\",1) or 1)\n    intercept=float(getattr(ds,\"RescaleIntercept\",0) or 0)\n    img=img*slope+intercept\n    lo,hi=np.percentile(img,[1,99])\n    img=np.clip((img-lo)/max(hi-lo,1e-6),0,1)\n    return img,ds\n\nfor study in studies:\n    study_dir=TRAIN_DIR/study\n    series_dirs=[p for p in study_dir.iterdir() if p.is_dir()]\n    series_dirs=series_dirs[:N_SERIES]\n    fig,axes=plt.subplots(1,len(series_dirs),figsize=(5*len(series_dirs),5))\n    if len(series_dirs)==1:axes=[axes]\n    for ax,sdir in zip(axes,series_dirs):\n        img,ds=read_middle_slice(sdir)\n        if img is None:\n            ax.axis(\"off\")\n            continue\n        desc=str(getattr(ds,\"SeriesDescription\",\"\"))\n        seq=str(getattr(ds,\"SequenceName\",\"\"))\n        spacing=getattr(ds,\"PixelSpacing\",\"\")\n        ax.imshow(img,cmap=\"gray\")\n        ax.set_title(f\"{study[:10]}...\\n{desc}\\n{seq}\\nPixelSpacing={spacing}\",fontsize=9)\n        ax.axis(\"off\")\n    plt.tight_layout()\n    plt.show()","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-08-16T02:20:37.7221Z","iopub.execute_input":"2026-08-16T02:20:37.722358Z","iopub.status.idle":"2026-08-16T02:20:47.996697Z","shell.execute_reply.started":"2026-08-16T02:20:37.722332Z","shell.execute_reply":"2026-08-16T02:20:47.99539Z"}},"outputs":[],"execution_count":null},{"cell_type":"code","source":"import os,re,random,warnings,hashlib,math\nfrom pathlib import Path\nfrom collections import Counter\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nimport pydicom\nwarnings.filterwarnings('ignore')\nSEED=1130\nrandom.seed(SEED);np.random.seed(SEED)\nTARGETS=['ACL','MCL','Medial Meniscus','Lateral Meniscus','Medial OA','Lateral OA','PF OA','Effusion','Synovitis',\"Baker's\",'Contusion','Fracture']\nHEADER_SAMPLE=5000\nGEOM_SAMPLE=300\nIMG_SAMPLE=500\nCROP_MM=130.0\nMODEL_IMG=336\n\ndef find_root():\n    for p in [Path('/kaggle/input/competitions/rsna-knee-abnormality-detection'),Path('/kaggle/input/rsna-knee-abnormality-detection')]:\n        if (p/'train.csv').exists() and (p/'train_series').exists():return p\n    for p in Path('/kaggle/input').glob('**/train.csv'):\n        if (p.parent/'train_series').exists():return p.parent\n    raise FileNotFoundError('Competition data not found')\n\ndef wilson(k,n,z=1.96):\n    if n<=0:return np.nan,np.nan\n    p=k/n;d=1+z*z/n;c=(p+z*z/(2*n))/d;h=z*np.sqrt(p*(1-p)/n+z*z/(4*n*n))/d\n    return max(0,c-h),min(1,c+h)\n\ndef normalize_report(s):\n    if not isinstance(s,str):return ''\n    s=s.lower().strip();s=re.sub(r'\\s+',' ',s);s=re.sub(r'[^\\w\\s]',' ',s);s=re.sub(r'\\s+',' ',s).strip();return s\n\ndef guess_language(text):\n    if not isinstance(text,str) or not text.strip():return '?'\n    t=text.lower()\n    if re.search(r'[\\u0370-\\u03ff]',t):return 'el'\n    if re.search(r'[\\u0400-\\u04ff]',t):return 'bg/ru'\n    words=re.findall(r'[a-zÀ-ÿıİ]+',t)\n    sets={'en':{'the','and','with','without','there','normal','no','knee','meniscus'},'es':{'del','los','las','con','sin','rodilla','menisco','rotura','ligamento'},'tr':{'ve','ile','diz','menisk','eklem','yok','bag','sivi','izlenmedi'},'hr':{'se','te','uz','bez','koljena','meniska','ligament','zgloba'},'de':{'der','die','das','mit','ohne','knie','meniskus','gelenk','band'},'nl':{'de','het','met','zonder','knie','meniscus','gewricht'},'fr':{'des','les','avec','sans','genou','menisque','ligament'}}\n    scores={k:sum(w in s for w in words) for k,s in sets.items()};best=max(scores,key=scores.get)\n    return best if scores[best]>=2 else '?'\n\ndef rank_corr(x,y):\n    x=np.asarray(x,float);y=np.asarray(y,float);ok=np.isfinite(x)&np.isfinite(y)\n    if ok.sum()<3:return np.nan\n    rx=pd.Series(x[ok]).rank(method='average').to_numpy();ry=pd.Series(y[ok]).rank(method='average').to_numpy()\n    if np.std(rx)==0 or np.std(ry)==0:return np.nan\n    return float(np.corrcoef(rx,ry)[0,1])\n\ndef safe_float(x):\n    try:return float(x)\n    except:return np.nan\n\ndef safe_first(seq):\n    try:return seq[0]\n    except:return np.nan\n\nROOT=find_root();TRAIN_DIR=ROOT/'train_series'\ntrain=pd.read_csv(ROOT/'train.csv',dtype={'StudyInstanceUID':str})\nseries=pd.read_csv(ROOT/'train_series.csv',dtype={'StudyInstanceUID':str,'SeriesInstanceUID':str})\nseries['Anatomical_Plane']=series['Anatomical_Plane'].astype(str)\nseries['Fluid_Sensitive']=pd.to_numeric(series['Fluid_Sensitive'],errors='coerce')\nseries['Fat_Suppression']=pd.to_numeric(series['Fat_Suppression'],errors='coerce')\ntrain['n_labels']=train[TARGETS].notna().sum(axis=1)\ntrain['n_positive_known']=train[TARGETS].fillna(0).sum(axis=1)\ntrain['report_chars']=train['Report'].fillna('').str.len()\ntrain['report_words']=train['Report'].fillna('').str.split().str.len()\ntrain['language']=train['Report'].fillna('').map(guess_language)\ntrain['report_norm']=train['Report'].fillna('').map(normalize_report)\ntrain['report_hash']=train['report_norm'].map(lambda s:hashlib.md5(s.encode()).hexdigest() if s else '')\nper_study=series.groupby('StudyInstanceUID').size().rename('n_series')\ntrain=train.merge(per_study,left_on='StudyInstanceUID',right_index=True,how='left')\ncomplete=train['n_labels'].eq(len(TARGETS));gold=train.loc[complete].copy()\nplt.rcParams.update({'figure.dpi':120,'font.size':9,'axes.spines.top':False,'axes.spines.right':False})\nprint('ROOT:',ROOT);print('train:',train.shape);print('train_series:',series.shape);print('studies:',train.StudyInstanceUID.nunique());print('series:',series.SeriesInstanceUID.nunique());print('fully annotated studies:',int(complete.sum()));print('reports:',int(train.Report.notna().sum()));print('missing reports:',int(train.Report.isna().sum()))\n\nann=train[TARGETS].notna().sum().sort_values();prev=[];lo=[];hi=[]\nfor t in ann.index:\n    y=train.loc[train[t].notna(),t].astype(int);prev.append(y.mean() if len(y) else np.nan);a,b=wilson(int(y.sum()),len(y));lo.append(a);hi.append(b)\nprev=pd.Series(prev,index=ann.index);lo=pd.Series(lo,index=ann.index);hi=pd.Series(hi,index=ann.index)\nfig,ax=plt.subplots(1,3,figsize=(19,4.8))\nax[0].barh(ann.index,ann.values);ax[0].set_xlabel('Annotated studies');ax[0].set_title('Label availability by target')\nfor i,v in enumerate(ann.values):ax[0].text(v+max(1,ann.max()*.01),i,str(v),va='center',fontsize=8)\ny=np.arange(len(prev));ax[1].errorbar(prev.values,y,xerr=[prev.values-lo.values,hi.values-prev.values],fmt='o',capsize=2);ax[1].set_yticks(y);ax[1].set_yticklabels(prev.index);ax[1].set_xlim(0,1);ax[1].set_xlabel('Positive rate among annotated');ax[1].set_title('Positive prevalence + 95% Wilson CI')\nvc=train.n_labels.value_counts().sort_index();ax[2].bar(vc.index.astype(str),vc.values);ax[2].set_xlabel('# available labels / study');ax[2].set_ylabel('Studies');ax[2].set_title('Supervision sparsity structure')\nplt.tight_layout();plt.show()\n\nfig,ax=plt.subplots(1,3,figsize=(19,4.6))\nnon=train.loc[~complete]\nax[0].boxplot([non.report_chars.values,gold.report_chars.values],labels=['Other','Fully annotated'],showfliers=False);ax[0].set_ylabel('Report characters');ax[0].set_title('Is the 58-study Gold subset typical?')\nax[1].boxplot([non.n_series.dropna().values,gold.n_series.dropna().values],labels=['Other','Fully annotated'],showfliers=False);ax[1].set_ylabel('Series / study');ax[1].set_title('Acquisition burden: Gold vs others')\nlangs=train.language.value_counts().head(8).index;allp=train.language.value_counts(normalize=True).reindex(langs).fillna(0);gp=gold.language.value_counts(normalize=True).reindex(langs).fillna(0);x=np.arange(len(langs));w=.38;ax[2].bar(x-w/2,allp.values,w,label='All');ax[2].bar(x+w/2,gp.values,w,label='Gold 58');ax[2].set_xticks(x);ax[2].set_xticklabels(langs,rotation=45);ax[2].set_ylabel('Fraction');ax[2].set_title('Language composition');ax[2].legend(fontsize=8)\nplt.tight_layout();plt.show()\n\nif len(gold)>=2:\n    corr=gold[TARGETS].astype(float).corr();poscount=gold[TARGETS].sum(axis=1)\n    pairs=[]\n    for i,a in enumerate(TARGETS):\n        for b in TARGETS[i+1:]:\n            ya=gold[a].astype(int);yb=gold[b].astype(int);joint=int(((ya==1)&(yb==1)).sum());den=int(((ya==1)|(yb==1)).sum());jac=joint/den if den else 0;pairs.append((f'{a} + {b}',joint,jac))\n    pairdf=pd.DataFrame(pairs,columns=['pair','joint','jaccard']).sort_values(['joint','jaccard'],ascending=False).head(12)\n    fig,ax=plt.subplots(1,3,figsize=(20,5.2))\n    im=ax[0].imshow(corr.values,vmin=-1,vmax=1,aspect='auto');ax[0].set_xticks(range(len(TARGETS)));ax[0].set_xticklabels(TARGETS,rotation=60,ha='right',fontsize=7);ax[0].set_yticks(range(len(TARGETS)));ax[0].set_yticklabels(TARGETS,fontsize=7);ax[0].set_title('Binary label correlation (Gold 58)');ax[0].grid(False);fig.colorbar(im,ax=ax[0],fraction=.046,pad=.04)\n    bins=np.arange(-.5,len(TARGETS)+1.5,1);ax[1].hist(poscount,bins=bins);ax[1].set_xlabel('# positive findings / study');ax[1].set_ylabel('Studies');ax[1].set_title('Multi-label disease burden')\n    z=pairdf.sort_values(['joint','jaccard']);ax[2].barh(z.pair,z.joint);ax[2].set_xlabel('Joint-positive studies');ax[2].set_title('Most common co-occurring findings')\n    plt.tight_layout();plt.show()\n\nlang=train.language.value_counts();dup=train.loc[train.report_hash!=''].groupby('report_hash').size().sort_values(ascending=False);dup_studies=int(dup[dup>1].sum()) if len(dup) else 0;dup_groups=int((dup>1).sum()) if len(dup) else 0\nfig,ax=plt.subplots(1,3,figsize=(19,4.5))\nax[0].bar(lang.index,lang.values);ax[0].tick_params(axis='x',rotation=45);ax[0].set_ylabel('Studies');ax[0].set_title('Estimated report language')\nvals=train.loc[train.report_chars>0,'report_chars'];ax[1].hist(vals,bins=45);ax[1].axvline(vals.median(),ls='--',lw=1);ax[1].set_xlabel('Characters');ax[1].set_ylabel('Studies');ax[1].set_title(f'Report length | median={vals.median():.0f}')\nif len(dup):\n    dv=dup.value_counts().sort_index();ax[2].bar(dv.index.astype(str),dv.values);ax[2].set_xlabel('Exact normalized report group size');ax[2].set_ylabel('Groups')\nax[2].set_title(f'Exact report duplicates\\n{dup_groups} duplicate groups / {dup_studies} studies')\nplt.tight_layout();plt.show()\nprint('\\nTop duplicated normalized reports:');print(dup.head(15).to_string())\n\nfig,ax=plt.subplots(1,3,figsize=(19,4.5))\nplane=series.Anatomical_Plane.value_counts();ax[0].bar(plane.index,plane.values);ax[0].set_ylabel('Series');ax[0].set_title('Anatomical plane')\ncross=pd.crosstab(series.Fluid_Sensitive,series.Fat_Suppression);x=np.arange(len(cross.index));bottom=np.zeros(len(cross))\nfor c in cross.columns:ax[1].bar(x,cross[c].values,bottom=bottom,label=f'FS={int(c)}');bottom+=cross[c].values\nax[1].set_xticks(x);ax[1].set_xticklabels([f'Fluid={int(v)}' for v in cross.index]);ax[1].set_ylabel('Series');ax[1].set_title('Fluid-sensitive × Fat suppression');ax[1].legend(fontsize=8)\nax[2].hist(per_study,bins=np.arange(per_study.min(),per_study.max()+2)-.5);ax[2].axvline(per_study.median(),ls='--',lw=1);ax[2].set_xlabel('Series / study');ax[2].set_ylabel('Studies');ax[2].set_title(f'Series count | median={per_study.median():.0f}')\nplt.tight_layout();plt.show()\n\npublic_slots=[('Sag Fluid','Sagittal',1),('Cor Fluid','Coronal',1),('Ax Fluid','Axial',1),('Sag Struct','Sagittal',0),('Cor Struct','Coronal',0),('Ax Struct','Axial',0)]\nslot_table=pd.DataFrame(index=train.StudyInstanceUID.astype(str).unique())\nfor name,pl,fluid in public_slots:\n    ids=series.loc[(series.Anatomical_Plane.eq(pl))&(series.Fluid_Sensitive.eq(fluid)),'StudyInstanceUID'].astype(str).unique();slot_table[name]=slot_table.index.isin(ids)\nslot_n=slot_table.sum(axis=1);slot_missing=(1-slot_table.mean()).sort_values(ascending=False);pattern=slot_table.astype(int).astype(str).agg(''.join,axis=1).value_counts().head(12)\nfig,ax=plt.subplots(1,3,figsize=(20,4.8))\nax[0].barh(slot_missing.index,slot_missing.values*100);ax[0].set_xlabel('Studies missing slot [%]');ax[0].set_title('Public 6-slot missingness')\nvc=slot_n.value_counts().sort_index();ax[1].bar(vc.index.astype(str),vc.values);ax[1].set_xlabel('# of 6 slots available');ax[1].set_ylabel('Studies');ax[1].set_title('Per-study protocol completeness')\nax[2].barh(pattern.index[::-1],pattern.values[::-1]);ax[2].set_xlabel('Studies');ax[2].set_title('Most common slot-presence patterns\\norder: SagF CorF AxF SagS CorS AxS')\nplt.tight_layout();plt.show()\n\nsample_series=series.sample(min(HEADER_SAMPLE,len(series)),random_state=SEED).reset_index(drop=True)\ndef inspect_series(row):\n    d=TRAIN_DIR/str(row.StudyInstanceUID)/str(row.SeriesInstanceUID)\n    try:files=[p for p in d.iterdir() if p.suffix.lower()=='.dcm']\n    except:return None\n    if not files:return None\n    f=files[len(files)//2]\n    tags=['PixelSpacing','Rows','Columns','SliceThickness','SpacingBetweenSlices','MagneticFieldStrength','Manufacturer','ManufacturerModelName','Laterality','ImageLaterality','SeriesDescription','RepetitionTime','EchoTime','ImagePositionPatient','ImageOrientationPatient','InstanceNumber']\n    try:\n        ds=pydicom.dcmread(f,stop_before_pixels=True,force=True,specific_tags=tags);ps=getattr(ds,'PixelSpacing',None);py=safe_float(ps[0]) if ps is not None and len(ps)>=2 else np.nan;px=safe_float(ps[1]) if ps is not None and len(ps)>=2 else np.nan;rows=safe_float(getattr(ds,'Rows',np.nan));cols=safe_float(getattr(ds,'Columns',np.nan));lat=str(getattr(ds,'Laterality','') or getattr(ds,'ImageLaterality','') or '')\n        fy=rows*py if np.isfinite(rows) and np.isfinite(py) else np.nan;fx=cols*px if np.isfinite(cols) and np.isfinite(px) else np.nan\n        return {'StudyInstanceUID':str(row.StudyInstanceUID),'SeriesInstanceUID':str(row.SeriesInstanceUID),'plane':row.Anatomical_Plane,'fluid':row.Fluid_Sensitive,'fat_suppression':row.Fat_Suppression,'n_slices':len(files),'pixel_y':py,'pixel_x':px,'rows':rows,'cols':cols,'fov_y_mm':fy,'fov_x_mm':fx,'fov_min_mm':min(fy,fx) if np.isfinite(fy) and np.isfinite(fx) else np.nan,'slice_thickness':safe_float(getattr(ds,'SliceThickness',np.nan)),'spacing_between':safe_float(getattr(ds,'SpacingBetweenSlices',np.nan)),'field_strength':safe_float(getattr(ds,'MagneticFieldStrength',np.nan)),'manufacturer':str(getattr(ds,'Manufacturer','Unknown')),'model':str(getattr(ds,'ManufacturerModelName','Unknown')),'laterality':lat,'description':str(getattr(ds,'SeriesDescription','')),'TR':safe_float(getattr(ds,'RepetitionTime',np.nan)),'TE':safe_float(getattr(ds,'EchoTime',np.nan))}\n    except:return None\nmeta=pd.DataFrame([x for x in (inspect_series(r) for _,r in sample_series.iterrows()) if x is not None]);meta.to_csv('/kaggle/working/eda_dicom_meta.csv',index=False)\nprint('\\nDICOM header sample:',len(meta),'series')\nfig,ax=plt.subplots(1,3,figsize=(19,4.5))\nq=meta.n_slices.quantile(.99);ax[0].hist(meta.loc[meta.n_slices<=q,'n_slices'],bins=40);ax[0].axvline(meta.n_slices.median(),ls='--',lw=1);ax[0].set_xlabel('Slices / series');ax[0].set_ylabel('Series');ax[0].set_title(f'Slice count | median={meta.n_slices.median():.0f}')\np=meta.pixel_x.dropna();q=p.quantile(.99);ax[1].hist(p[p<=q],bins=40);ax[1].axvline(p.median(),ls='--',lw=1);ax[1].axvline(CROP_MM/MODEL_IMG,ls=':',lw=1,label=f'{CROP_MM/MODEL_IMG:.3f} mm/px @ {MODEL_IMG}px');ax[1].set_xlabel('Pixel spacing [mm/pixel]');ax[1].set_ylabel('Series');ax[1].set_title('Native in-plane resolution');ax[1].legend(fontsize=8)\nf=meta.fov_min_mm.dropna();loq,hiq=f.quantile([.01,.99]);ax[2].hist(f[(f>=loq)&(f<=hiq)],bins=40);ax[2].axvline(CROP_MM,ls='--',lw=1,label='130 mm crop');ax[2].set_xlabel('Minimum physical FOV [mm]');ax[2].set_ylabel('Series');ax[2].set_title('Physical field of view');ax[2].legend(fontsize=8)\nplt.tight_layout();plt.show()\n\nfig,ax=plt.subplots(1,3,figsize=(19,4.6))\nfor pl in ['Sagittal','Coronal','Axial']:\n    x=meta.loc[meta.plane.eq(pl),'pixel_x'].dropna().values\n    if len(x):ax[0].boxplot(x,positions=[['Sagittal','Coronal','Axial'].index(pl)+1],widths=.55,showfliers=False)\nax[0].set_xticks([1,2,3]);ax[0].set_xticklabels(['Sagittal','Coronal','Axial']);ax[0].set_ylabel('Pixel spacing [mm]');ax[0].set_title('In-plane resolution by plane')\nfor pl in ['Sagittal','Coronal','Axial']:\n    x=meta.loc[meta.plane.eq(pl),'slice_thickness'].dropna().values\n    if len(x):ax[1].boxplot(x,positions=[['Sagittal','Coronal','Axial'].index(pl)+1],widths=.55,showfliers=False)\nax[1].set_xticks([1,2,3]);ax[1].set_xticklabels(['Sagittal','Coronal','Axial']);ax[1].set_ylabel('Slice thickness [mm]');ax[1].set_title('Through-plane resolution by plane')\nfor pl in ['Sagittal','Coronal','Axial']:\n    x=meta.loc[meta.plane.eq(pl),'n_slices'].dropna().values\n    if len(x):ax[2].boxplot(x,positions=[['Sagittal','Coronal','Axial'].index(pl)+1],widths=.55,showfliers=False)\nax[2].set_xticks([1,2,3]);ax[2].set_xticklabels(['Sagittal','Coronal','Axial']);ax[2].set_ylabel('Slices / series');ax[2].set_title('Volume depth by plane')\nplt.tight_layout();plt.show()\n\nfig,ax=plt.subplots(1,3,figsize=(20,4.6))\nm=meta.manufacturer.replace({'':'Unknown'}).value_counts().head(10);ax[0].barh(m.index[::-1],m.values[::-1]);ax[0].set_xlabel('Series');ax[0].set_title('Scanner manufacturer')\nfs=meta.field_strength.dropna().round(1).value_counts().sort_index();ax[1].bar(fs.index.astype(str),fs.values);ax[1].set_xlabel('Tesla');ax[1].set_ylabel('Series');ax[1].set_title('Magnetic field strength')\nmodel=meta.model.replace({'':'Unknown'}).value_counts().head(12);ax[2].barh(model.index[::-1],model.values[::-1]);ax[2].set_xlabel('Series');ax[2].set_title('Scanner model')\nplt.tight_layout();plt.show()\n\ngeom_series=series.sample(min(GEOM_SAMPLE,len(series)),random_state=SEED+1).reset_index(drop=True)\ndef geometry_qc(row):\n    d=TRAIN_DIR/str(row.StudyInstanceUID)/str(row.SeriesInstanceUID)\n    files=sorted(d.glob('*.dcm'))\n    if len(files)<3:return None\n    rec=[]\n    for j,f in enumerate(files):\n        try:\n            ds=pydicom.dcmread(f,stop_before_pixels=True,force=True,specific_tags=['ImagePositionPatient','ImageOrientationPatient','InstanceNumber']);ipp=getattr(ds,'ImagePositionPatient',None);iop=getattr(ds,'ImageOrientationPatient',None);inst=safe_float(getattr(ds,'InstanceNumber',np.nan))\n            if ipp is not None and iop is not None and len(ipp)>=3 and len(iop)>=6:\n                ipp=np.asarray(ipp,float);iop=np.asarray(iop,float);pos=float(np.dot(ipp,np.cross(iop[:3],iop[3:6])))\n            else:pos=np.nan\n            rec.append((j,inst,pos))\n        except:rec.append((j,np.nan,np.nan))\n    a=np.asarray(rec,float);valid=np.isfinite(a[:,2]);coverage=float(valid.mean());fname_corr=rank_corr(a[:,0],a[:,2]);inst_corr=rank_corr(a[:,1],a[:,2]);sp=np.diff(np.sort(a[valid,2])) if valid.sum()>=3 else np.array([]);sp=np.abs(sp[np.abs(sp)>1e-6]);spacing_med=float(np.median(sp)) if len(sp) else np.nan;spacing_cv=float(np.std(sp)/np.mean(sp)) if len(sp)>1 and np.mean(sp)>0 else np.nan\n    return {'StudyInstanceUID':str(row.StudyInstanceUID),'SeriesInstanceUID':str(row.SeriesInstanceUID),'plane':row.Anatomical_Plane,'n':len(files),'geometry_coverage':coverage,'filename_rho':fname_corr,'instance_rho':inst_corr,'slice_spacing_mm':spacing_med,'spacing_cv':spacing_cv}\ngeom=pd.DataFrame([x for x in (geometry_qc(r) for _,r in geom_series.iterrows()) if x is not None]);geom.to_csv('/kaggle/working/eda_geometry_qc.csv',index=False)\nfig,ax=plt.subplots(1,3,figsize=(19,4.5))\nax[0].hist(geom.geometry_coverage*100,bins=20);ax[0].set_xlabel('Slices with IOP+IPP [%]');ax[0].set_ylabel('Series');ax[0].set_title('Can physical ordering be recovered?')\nfc=geom.filename_rho.abs().dropna();ic=geom.instance_rho.abs().dropna();ax[1].boxplot([fc.values,ic.values],labels=['Filename order','InstanceNumber']);ax[1].set_ylabel('|Spearman ρ| vs physical position');ax[1].set_ylim(0,1.05);ax[1].set_title('Which ordering signal matches anatomy?')\ncv=geom.spacing_cv.replace([np.inf,-np.inf],np.nan).dropna();q=cv.quantile(.99) if len(cv) else 1\nif len(cv):ax[2].hist(cv[cv<=q],bins=30)\nax[2].set_xlabel('CV of adjacent physical slice spacing');ax[2].set_ylabel('Series');ax[2].set_title('Slice-spacing regularity')\nplt.tight_layout();plt.show()\n\nlat=meta.laterality.fillna('').str.upper().str.strip();tagged=lat.isin(['L','R','LEFT','RIGHT']);aspect=meta.cols/meta.rows;fov_ratio=meta.fov_x_mm/meta.fov_y_mm\nfig,ax=plt.subplots(1,3,figsize=(19,4.5))\nax[0].bar(['Tagged','Missing/other'],[int(tagged.sum()),int((~tagged).sum())]);ax[0].set_ylabel('Series');ax[0].set_title('Laterality tag coverage')\nar=aspect.replace([np.inf,-np.inf],np.nan).dropna();ax[1].hist(ar,bins=35);ax[1].axvline(1,ls='--',lw=1);ax[1].set_xlabel('Columns / Rows');ax[1].set_ylabel('Series');ax[1].set_title('Matrix aspect ratio')\nfr=fov_ratio.replace([np.inf,-np.inf],np.nan).dropna();ax[2].hist(fr,bins=35);ax[2].axvline(1,ls='--',lw=1);ax[2].set_xlabel('FOVx / FOVy');ax[2].set_ylabel('Series');ax[2].set_title('Physical FOV aspect ratio')\nplt.tight_layout();plt.show()\n\nimg_rows=[]\nfor _,row in sample_series.sample(min(IMG_SAMPLE,len(sample_series)),random_state=SEED+2).iterrows():\n    d=TRAIN_DIR/str(row.StudyInstanceUID)/str(row.SeriesInstanceUID);files=sorted(d.glob('*.dcm'))\n    if not files:continue\n    try:\n        ds=pydicom.dcmread(files[len(files)//2],force=True);a=ds.pixel_array.astype(np.float32);a=a*safe_float(getattr(ds,'RescaleSlope',1))+safe_float(getattr(ds,'RescaleIntercept',0));p1,p50,p99=np.percentile(a,[1,50,99]);img_rows.append({'plane':row.Anatomical_Plane,'p1':p1,'p50':p50,'p99':p99,'dynamic':p99-p1,'median':p50})\n    except:pass\nimgmeta=pd.DataFrame(img_rows)\nfig,ax=plt.subplots(1,3,figsize=(19,4.5))\nif len(imgmeta):\n    dyn=np.log1p(imgmeta.dynamic.clip(lower=0));ax[0].hist(dyn,bins=40);ax[0].set_xlabel('log(1 + p99-p1)');ax[0].set_ylabel('Slices');ax[0].set_title('Raw MRI intensity scale varies strongly')\n    planes=['Sagittal','Coronal','Axial'];vals=[np.log1p(imgmeta.loc[imgmeta.plane.eq(p),'dynamic'].clip(lower=0).values) for p in planes];ax[1].boxplot(vals,labels=planes,showfliers=False);ax[1].set_ylabel('log dynamic range');ax[1].set_title('Intensity scale by plane')\n    med=np.log1p(np.abs(imgmeta['median']));ax[2].hist(med,bins=40);ax[2].set_xlabel('log(1 + |median intensity|)');ax[2].set_ylabel('Slices');ax[2].set_title('Raw intensity offset variability')\nplt.tight_layout();plt.show()\n\ndef physical_order(series_dir):\n    files=list(series_dir.glob('*.dcm'));rows=[]\n    for f in files:\n        try:\n            ds=pydicom.dcmread(f,stop_before_pixels=True,force=True,specific_tags=['ImagePositionPatient','ImageOrientationPatient','InstanceNumber']);ipp=getattr(ds,'ImagePositionPatient',None);iop=getattr(ds,'ImageOrientationPatient',None)\n            if ipp is not None and iop is not None:\n                ipp=np.asarray(ipp,float);iop=np.asarray(iop,float);key=float(np.dot(ipp,np.cross(iop[:3],iop[3:6])))\n            else:key=safe_float(getattr(ds,'InstanceNumber',np.nan))\n        except:key=np.nan\n        rows.append((key,f))\n    if rows and all(np.isfinite(k) for k,_ in rows):return [f for k,f in sorted(rows,key=lambda x:x[0])]\n    return sorted(files)\n\ndef read_middle(series_dir):\n    files=physical_order(series_dir)\n    if not files:return None,None\n    for off in [0,-1,1,-2,2]:\n        i=int(np.clip(len(files)//2+off,0,len(files)-1))\n        try:\n            ds=pydicom.dcmread(files[i],force=True);a=ds.pixel_array.astype(np.float32);a=a*safe_float(getattr(ds,'RescaleSlope',1))+safe_float(getattr(ds,'RescaleIntercept',0));p1,p99=np.percentile(a,[1,99]);a=np.clip((a-p1)/max(p99-p1,1e-6),0,1);return a,ds\n        except:pass\n    return None,None\nthree=series.groupby('StudyInstanceUID').Anatomical_Plane.nunique();candidates=three[three>=3].index.astype(str).tolist();rng=np.random.default_rng(SEED);study=str(rng.choice(candidates));g=series[series.StudyInstanceUID.astype(str).eq(study)]\nfig,ax=plt.subplots(1,3,figsize=(15,5))\nfor j,pl in enumerate(['Sagittal','Coronal','Axial']):\n    cand=g[g.Anatomical_Plane.eq(pl)]\n    if cand.empty:ax[j].axis('off');continue\n    row=cand.sort_values(['Fluid_Sensitive','Fat_Suppression'],ascending=False).iloc[0];img,ds=read_middle(TRAIN_DIR/study/str(row.SeriesInstanceUID))\n    if img is None:ax[j].text(.5,.5,'Decode failed',ha='center',va='center');ax[j].axis('off');continue\n    ax[j].imshow(img,cmap='gray');ax[j].set_title(f\"{pl}\\nFluid={row.Fluid_Sensitive}  FS={row.Fat_Suppression}\\n{str(getattr(ds,'SeriesDescription',''))[:32]}\",fontsize=9);ax[j].axis('off')\nplt.suptitle(f'One study, three anatomical planes | {study}',fontsize=11);plt.tight_layout();plt.show()\n\nsummary=[]\ndef add(metric,value):summary.append((metric,value))\nadd('Studies',len(train));add('Series',len(series));add('Fully annotated studies',int(complete.sum()));add('Median labels available / study',float(train.n_labels.median()));add('Exact duplicate report groups',dup_groups);add('Studies inside duplicate report groups',dup_studies);add('Studies with all 6 public slots',int((slot_n==6).sum()));add('All-6-slot rate [%]',round(float((slot_n==6).mean()*100),1));add('Header sample series',len(meta));add('Laterality tag coverage [%]',round(float(tagged.mean()*100),1));add('Series with min FOV < 130 mm [%]',round(float((meta.fov_min_mm<CROP_MM).mean()*100),1));add('Median native pixel spacing [mm]',round(float(meta.pixel_x.median()),3));add('Model target spacing at 336 px [mm]',round(CROP_MM/MODEL_IMG,3));add('Geometry QC series',len(geom));add('Median geometry coverage [%]',round(float(geom.geometry_coverage.median()*100),1));add('Median |rho| filename→physical',round(float(geom.filename_rho.abs().median()),3));add('Median |rho| InstanceNumber→physical',round(float(geom.instance_rho.abs().median()),3));add('Median slice-spacing CV',round(float(geom.spacing_cv.median()),3))\nsummary=pd.DataFrame(summary,columns=['EDA metric','value']);print('\\n=== MODEL-RELEVANT EDA SUMMARY ===');print(summary.to_string(index=False));summary.to_csv('/kaggle/working/eda_summary.csv',index=False)\nprint('\\nSaved: /kaggle/working/eda_dicom_meta.csv, eda_geometry_qc.csv, eda_summary.csv')\n","metadata":{"trusted":true,"execution":{"iopub.status.busy":"2026-08-16T02:20:47.998982Z","iopub.execute_input":"2026-08-16T02:20:47.999276Z","iopub.status.idle":"2026-08-16T02:24:31.327938Z","shell.execute_reply.started":"2026-08-16T02:20:47.999247Z","shell.execute_reply":"2026-08-16T02:24:31.326798Z"}},"outputs":[],"execution_count":null},{"cell_type":"markdown","source":"Exploratory data analysis revealed substantial heterogeneity in both supervision and MRI acquisition. The dataset contains 4,407 studies and 24,371 MRI series, but only 58 studies are fully annotated for all 12 target abnormalities. The median number of available manual labels per study is therefore zero, confirming that conventional fully supervised learning is not feasible for the majority of the training set and motivating the use of report-derived weak supervision.\nThe analysis also identified 52 groups of duplicated normalized reports involving 223 studies. Because identical reports can yield nearly identical pseudo-labels, randomly distributing these studies across training and validation sets could lead to information leakage. This finding supports grouping studies with identical normalized reports within the same cross-validation fold.\nConsiderable variation was also observed in MRI acquisition protocols. Only 566 studies (12.8%) contained all six predefined combinations of anatomical plane and sequence type, indicating that complete multi-sequence coverage is uncommon. Thus, a model should not assume that all six MRI slots are available for every study. Instead, explicitly representing sequence availability and allowing missing slots is necessary for robust inference across heterogeneous acquisition protocols.\nThe DICOM metadata further support the adopted spatial preprocessing strategy. In a sample of 5,000 series, the median native in-plane pixel spacing was 0.312 mm/pixel, whereas resampling a 130-mm physical crop to 336 × 336 pixels corresponds to 0.387 mm/pixel. Importantly, only 0.4% of sampled series had a minimum physical field of view smaller than 130 mm, indicating that a 130-mm crop can be applied to nearly all acquisitions without exceeding the available field of view. These findings support physical-scale normalization before resizing rather than directly resizing the original pixel matrix.\nSlice-order analysis showed an especially strong justification for geometry-aware preprocessing. Across 300 sampled series, the median availability of ImagePositionPatient and ImageOrientationPatient information was 100%. File-name order exhibited only a weak relationship with anatomical position, with a median absolute Spearman correlation of 0.122, whereas InstanceNumber showed a median absolute correlation of 1.000 with physical slice position. In addition, the median coefficient of variation of physical slice spacing was 0.000, indicating highly regular spacing for most evaluated series. These results demonstrate that file names should not be used to determine slice order. Physical DICOM geometry provides the most rigorous ordering criterion, while InstanceNumber represents a highly reliable fallback in this dataset.\nLaterality metadata were available for only 50.0% of sampled series. Therefore, relying exclusively on the DICOM laterality tag would leave approximately half of the acquisitions unresolved. This is particularly relevant for medial–lateral targets such as the medial and lateral menisci and tibiofemoral osteoarthritis, and supports supplementing laterality tags with patient-coordinate geometry.\nOverall, the EDA indicates that the principal challenges of this dataset are sparse supervision, duplicated report-derived labels, incomplete MRI sequence coverage, and heterogeneous imaging geometry. These observations directly motivate leakage-aware grouped cross-validation, missing-slot-aware multi-view modeling, physical-scale normalization, geometry-based slice ordering, and geometry-assisted laterality normalization.","metadata":{}},{"cell_type":"code","source":"from __future__ import annotations\n\nimport json\nfrom pathlib import Path\nfrom itertools import combinations\n\nimport numpy as np\nimport pandas as pd\nfrom sklearn.metrics import roc_auc_score\n\nSEED = 1130\n\nTARGETS = [\n    \"ACL\",\"MCL\",\"Medial Meniscus\",\"Lateral Meniscus\",\n    \"Medial OA\",\"Lateral OA\",\"PF OA\",\"Effusion\",\n    \"Synovitis\",\"Baker's\",\"Contusion\",\"Fracture\"\n]\nARMS = [\"dino\",\"a5\",\"rad\",\"b3\"]\n\nGRID_STEP = 0.05\nDINO_MIN = 0.50\nA5_MAX = 0.30\nRAD_MAX = 0.40\nB3_MAX = 0.20\n\nGOLD_SCORE_WEIGHT = 0.50\nPSEUDO_SCORE_WEIGHT = 0.35\nFOLD_NONLOSS_WEIGHT = 0.10\nFOLD_GAIN_WEIGHT = 0.05\nDINO_PRIOR_PENALTY = 0.005\n\nTARGET_SHRINK_TO_GLOBAL = 0.35\nFOLD_NONLOSS_TOL = 0.002\nFOLD_GAIN_SCALE = 0.01\nGOLD_BOOTSTRAPS = 500\n\nOUT = Path(\"/kaggle/working/rsna_oof_blend_optimizer\")\nOUT.mkdir(parents=True, exist_ok=True)\n\ndef find_root():\n    cands = [\n        Path(\"/kaggle/input/competitions/rsna-knee-abnormality-detection\"),\n        Path(\"/kaggle/input/rsna-knee-abnormality-detection\"),\n    ]\n    for p in cands:\n        if (p/\"train.csv\").is_file():\n            return p\n    for p in Path(\"/kaggle/input\").glob(\"**/train.csv\"):\n        if (p.parent/\"train_series.csv\").is_file():\n            return p.parent\n    raise FileNotFoundError(\"RSNA competition root not found\")\n\ndef _all_named(filename):\n    out = []\n    for root in (Path(\"/kaggle/input\"), Path(\"/kaggle/working\")):\n        if root.exists():\n            out.extend(root.glob(f\"**/{filename}\"))\n    return sorted(set(out))\n\ndef _near_text(p):\n    parts = [p.parent.name.lower()]\n    for name in (\"training_config.json\", \"rad_heads_manifest.json\", \"manifest.json\"):\n        q = p.parent / name\n        if q.is_file():\n            try:\n                parts.append(q.read_text(errors=\"ignore\").lower())\n            except Exception:\n                pass\n    return \" \".join(parts)\n\ndef choose_file(filename, include_keywords, exclude_keywords=()):\n    files = _all_named(filename)\n    scored = []\n    for p in files:\n        near = _near_text(p)\n        full = str(p).lower()\n\n        # Arm identity is judged primarily from the immediate output folder/config.\n        # Dataset ancestor names such as \"...a5-b3-rad...\" must not exclude a valid Rad OOF.\n        score = 0.0\n        score += sum(100.0 for x in include_keywords if x.lower() in near)\n        score += sum(10.0 for x in include_keywords if x.lower() in full)\n        score -= sum(120.0 for x in exclude_keywords if x.lower() in near)\n        score -= len(p.parts) * 0.001\n        scored.append((score, p, near[:300]))\n\n    if not scored:\n        raise FileNotFoundError(f\"{filename} was not found anywhere under /kaggle/input or /kaggle/working\")\n\n    scored.sort(key=lambda x: (x[0], str(x[1])), reverse=True)\n    best = scored[0]\n\n    if not any(x.lower() in _near_text(best[1]) or x.lower() in str(best[1]).lower()\n               for x in include_keywords):\n        print(f\"WARNING: weak match for {include_keywords}: {best[1]}\")\n\n    return best[1]\n\ndef _show_oof_candidates():\n    rows = []\n    for fn in (\"oof_rankmean.csv\", \"oof_predictions.csv\"):\n        for p in _all_named(fn):\n            rows.append({\n                \"file\": fn,\n                \"parent\": p.parent.name,\n                \"path\": str(p),\n                \"near\": _near_text(p)[:180],\n            })\n    if rows:\n        print(\"\\nOOF candidates visible to this notebook:\")\n        for r in rows:\n            print(f\"  {r['file']:<20} | {r['parent']:<34} | {r['path']}\")\n    else:\n        print(\"\\nNo OOF CSV candidates are visible under /kaggle/input or /kaggle/working.\")\n\ndef discover():\n    try:\n        dino = choose_file(\"oof_rankmean.csv\", [\"dinov2\", \"20member\"])\n        a5 = choose_file(\"oof_predictions.csv\", [\"a5\", \"convnext\"], [\"b3\", \"radimagenet\"])\n        rad = choose_file(\"oof_predictions.csv\", [\"radimagenet\"], [\"a5_convnext\", \"b3_mil\"])\n        b3 = choose_file(\"oof_predictions.csv\", [\"b3\", \"mil\"], [\"a5_convnext\", \"radimagenet\"])\n    except FileNotFoundError:\n        _show_oof_candidates()\n        raise\n\n    pseudo = dino.parent / \"pseudo_labels.csv\"\n    if not pseudo.is_file():\n        pseudo = choose_file(\"pseudo_labels.csv\", [\"dinov2\", \"20member\"])\n\n    return {\"dino\": dino, \"a5\": a5, \"rad\": rad, \"b3\": b3, \"pseudo\": pseudo}\n\ndef auc_binary(y, p):\n    y = np.asarray(y, dtype=float)\n    p = np.asarray(p, dtype=float)\n    ok = np.isfinite(y) & np.isfinite(p)\n    y, p = y[ok], p[ok]\n    if len(y) < 2:\n        return np.nan\n    y = (y >= 0.5).astype(np.int8)\n    if len(np.unique(y)) < 2:\n        return np.nan\n    return float(roc_auc_score(y, p))\n\ndef auc_soft(y, p):\n    \"\"\"\n    ROC-AUC for probabilistic/soft pseudo labels in [0,1].\n\n    Each sample contributes fractionally to the positive and negative classes:\n      positive weight = y * confidence\n      negative weight = (1-y) * confidence\n    where confidence = 2*abs(y-0.5).\n\n    Thus y=0.5 (fully ambiguous) contributes zero weight, while 0/1 labels\n    reduce to ordinary weighted binary ROC-AUC.\n    \"\"\"\n    y = np.asarray(y, dtype=float)\n    p = np.asarray(p, dtype=float)\n    ok = np.isfinite(y) & np.isfinite(p)\n    y, p = y[ok], p[ok]\n    if len(y) < 2:\n        return np.nan\n\n    y = np.clip(y, 0.0, 1.0)\n    conf = 2.0 * np.abs(y - 0.5)\n    keep = conf > 1e-8\n    y, p, conf = y[keep], p[keep], conf[keep]\n    if len(y) < 2:\n        return np.nan\n\n    labels = np.concatenate([\n        np.ones(len(y), dtype=np.int8),\n        np.zeros(len(y), dtype=np.int8),\n    ])\n    scores = np.concatenate([p, p])\n    weights = np.concatenate([\n        y * conf,\n        (1.0 - y) * conf,\n    ])\n\n    if weights[:len(y)].sum() <= 0 or weights[len(y):].sum() <= 0:\n        return np.nan\n\n    return float(roc_auc_score(labels, scores, sample_weight=weights))\n\ndef fold_rank(df):\n    out=df[[\"StudyInstanceUID\",\"fold\"]].copy()\n    for t in TARGETS:\n        out[t]=df.groupby(\"fold\",dropna=False)[t].rank(method=\"average\",pct=True)\n    return out\n\ndef make_grid():\n    units=int(round(1/GRID_STEP))\n    ws=[]\n    for d in range(units+1):\n        for a in range(units-d+1):\n            for r in range(units-d-a+1):\n                b=units-d-a-r\n                w=np.array([d,a,r,b],dtype=float)/units\n                if w[0] < DINO_MIN: continue\n                if w[1] > A5_MAX: continue\n                if w[2] > RAD_MAX: continue\n                if w[3] > B3_MAX: continue\n                ws.append(w)\n    W=np.unique(np.round(np.vstack(ws),8),axis=0)\n    return W\n\nROOT=find_root()\nFILES=discover()\nprint(\"OOF sources:\")\nfor k,v in FILES.items():\n    print(f\"  {k:>6}: {v}\")\n\nfor arm in (\"dino\", \"a5\", \"rad\", \"b3\"):\n    z = pd.read_csv(FILES[arm], nrows=3)\n    required = {\"StudyInstanceUID\", \"fold\", *TARGETS}\n    missing = required - set(z.columns)\n    if missing:\n        raise ValueError(f\"{arm}: selected OOF file has wrong schema: {FILES[arm]} | missing={sorted(missing)}\")\n\ntrain=pd.read_csv(ROOT/\"train.csv\",dtype={\"StudyInstanceUID\":str})\npseudo=pd.read_csv(FILES[\"pseudo\"],dtype={\"StudyInstanceUID\":str})\n\nframes={}\nfor arm in ARMS:\n    z=pd.read_csv(FILES[arm],dtype={\"StudyInstanceUID\":str})\n    need=[\"StudyInstanceUID\",\"fold\"]+TARGETS\n    miss=[c for c in need if c not in z.columns]\n    if miss:\n        raise ValueError(f\"{arm}: missing columns {miss}\")\n    z=z[need].copy()\n    z[\"StudyInstanceUID\"]=z[\"StudyInstanceUID\"].astype(str)\n    frames[arm]=fold_rank(z)\n\nbase=frames[\"dino\"][[\"StudyInstanceUID\",\"fold\"]].rename(columns={\"fold\":\"meta_fold\"}).copy()\nfor arm in ARMS:\n    z=frames[arm].rename(columns={t:f\"{arm}__{t}\" for t in TARGETS})\n    z=z.rename(columns={\"fold\":f\"{arm}__fold\"})\n    base=base.merge(z,on=\"StudyInstanceUID\",how=\"inner\",validate=\"one_to_one\")\n\ngold_mask=train[TARGETS].notna().all(axis=1)\ngold=train.loc[gold_mask,[\"StudyInstanceUID\"]+TARGETS].copy()\ngold[\"StudyInstanceUID\"]=gold[\"StudyInstanceUID\"].astype(str)\n\npcols=[\"StudyInstanceUID\"]+[t for t in TARGETS if t in pseudo.columns]\npseudo=pseudo[pcols].copy()\npseudo[\"StudyInstanceUID\"]=pseudo[\"StudyInstanceUID\"].astype(str)\n\ndata=base.merge(pseudo,on=\"StudyInstanceUID\",how=\"left\",validate=\"one_to_one\",suffixes=(\"\",\"__pseudo\"))\ngold_ren=gold.rename(columns={t:f\"gold__{t}\" for t in TARGETS})\ndata=data.merge(gold_ren,on=\"StudyInstanceUID\",how=\"left\",validate=\"one_to_one\")\ndata[\"is_gold\"]=data[\"StudyInstanceUID\"].isin(set(gold[\"StudyInstanceUID\"]))\n\nprint(\"\\nPseudo-label diagnostics:\")\nfor t in TARGETS:\n    v = pd.to_numeric(data[t], errors=\"coerce\").to_numpy(float)\n    v = v[np.isfinite(v)]\n    if len(v):\n        frac_soft = float(np.mean((v > 0.0) & (v < 1.0)))\n        print(f\"  {t:<18} n={len(v):4d}  min={v.min():.3f}  max={v.max():.3f}  soft_fraction={frac_soft:.1%}\")\n\nfold_agreement=[]\nfor arm in [\"a5\",\"rad\",\"b3\"]:\n    c=(data[\"meta_fold\"].to_numpy()==data[f\"{arm}__fold\"].to_numpy())\n    fold_agreement.append({\"arm\":arm,\"fold_agreement_with_dino\":float(np.mean(c))})\npd.DataFrame(fold_agreement).to_csv(OUT/\"fold_agreement.csv\",index=False)\nprint(\"\\nFold agreement with DINO:\")\nprint(pd.DataFrame(fold_agreement).to_string(index=False))\n\ncorr_rows=[]\nfor t in TARGETS:\n    for x,y in combinations(ARMS,2):\n        c=data[[f\"{x}__{t}\",f\"{y}__{t}\"]].corr(method=\"spearman\").iloc[0,1]\n        corr_rows.append({\"target\":t,\"arm1\":x,\"arm2\":y,\"spearman\":float(c)})\ncorr_df=pd.DataFrame(corr_rows)\ncorr_df.to_csv(OUT/\"arm_rank_correlations.csv\",index=False)\n\nW=make_grid()\nprint(f\"\\nWeight grid: {len(W)} candidates\")\nprint(f\"Constraints: DINO>={DINO_MIN:.2f}, A5<={A5_MAX:.2f}, Rad<={RAD_MAX:.2f}, B3<={B3_MAX:.2f}\")\n\nall_candidates=[]\nper_target={}\n\nfor t in TARGETS:\n    X=np.column_stack([data[f\"{a}__{t}\"].to_numpy(float) for a in ARMS])\n    pseudo_y=data[t].to_numpy(float)\n    gold_y=data[f\"gold__{t}\"].to_numpy(float)\n    non_gold=(~data[\"is_gold\"]).to_numpy()\n    gold_rows=data[\"is_gold\"].to_numpy()\n    folds=data[\"meta_fold\"].to_numpy()\n\n    dino=X[:,0]\n    dino_fold_auc={}\n    for f in sorted(pd.Series(folds).dropna().unique()):\n        m=(folds==f)&non_gold\n        dino_fold_auc[f]=auc_soft(pseudo_y[m],dino[m])\n\n    rows=[]\n    for ci,w in enumerate(W):\n        p=X@w\n        pa=auc_soft(pseudo_y[non_gold],p[non_gold])\n        ga=auc_binary(gold_y[gold_rows],p[gold_rows])\n\n        gains=[]\n        for f in sorted(pd.Series(folds).dropna().unique()):\n            m=(folds==f)&non_gold\n            ca=auc_soft(pseudo_y[m],p[m])\n            da=dino_fold_auc[f]\n            if np.isfinite(ca) and np.isfinite(da):\n                gains.append(ca-da)\n\n        gain_mean=float(np.mean(gains)) if gains else 0.0\n        nonloss=float(np.mean(np.asarray(gains)>=-FOLD_NONLOSS_TOL)) if gains else 0.0\n        gain_score=float(np.clip(0.5+gain_mean/(2*FOLD_GAIN_SCALE),0,1))\n        obj=(\n            GOLD_SCORE_WEIGHT*(ga if np.isfinite(ga) else 0.5)\n            +PSEUDO_SCORE_WEIGHT*(pa if np.isfinite(pa) else 0.5)\n            +FOLD_NONLOSS_WEIGHT*nonloss\n            +FOLD_GAIN_WEIGHT*gain_score\n            -DINO_PRIOR_PENALTY*(1-w[0])\n        )\n        row={\n            \"target\":t,\"candidate\":ci,\n            \"dino\":w[0],\"a5\":w[1],\"rad\":w[2],\"b3\":w[3],\n            \"gold_auc\":ga,\"pseudo_auc\":pa,\n            \"fold_gain_mean\":gain_mean,\"fold_nonloss_rate\":nonloss,\n            \"objective\":obj\n        }\n        rows.append(row)\n        all_candidates.append(row)\n\n    td=pd.DataFrame(rows).sort_values([\"objective\",\"gold_auc\",\"pseudo_auc\"],ascending=False).reset_index(drop=True)\n    per_target[t]=td\n\ncand_df=pd.DataFrame(all_candidates)\ncand_df.to_csv(OUT/\"all_target_candidates.csv\",index=False)\n\nglobal_df=(cand_df.groupby(\"candidate\",as_index=False)\n    .agg(\n        dino=(\"dino\",\"first\"),a5=(\"a5\",\"first\"),rad=(\"rad\",\"first\"),b3=(\"b3\",\"first\"),\n        macro_gold_auc=(\"gold_auc\",\"mean\"),\n        macro_pseudo_auc=(\"pseudo_auc\",\"mean\"),\n        mean_fold_gain=(\"fold_gain_mean\",\"mean\"),\n        mean_fold_nonloss=(\"fold_nonloss_rate\",\"mean\"),\n        macro_objective=(\"objective\",\"mean\")\n    )\n    .sort_values([\"macro_objective\",\"macro_gold_auc\",\"macro_pseudo_auc\"],ascending=False)\n    .reset_index(drop=True)\n)\nglobal_df.to_csv(OUT/\"global_candidate_ranking.csv\",index=False)\ng=global_df.iloc[0]\nglobal_w=np.array([g[a] for a in ARMS],dtype=float)\n\nprint(\"\\nSelected global weight:\")\nprint({a:round(float(global_w[i]),4) for i,a in enumerate(ARMS)})\nprint(f\"macro Gold={g.macro_gold_auc:.4f} | pseudo={g.macro_pseudo_auc:.4f} | objective={g.macro_objective:.4f}\")\n\nrng=np.random.default_rng(SEED)\nfinal_rows=[]\ntop_rows=[]\n\nfor t in TARGETS:\n    td=per_target[t]\n    top_rows.append(td.head(10))\n    raw=td.iloc[0]\n    raw_w=np.array([raw[a] for a in ARMS],dtype=float)\n    final_w=(1-TARGET_SHRINK_TO_GLOBAL)*raw_w + TARGET_SHRINK_TO_GLOBAL*global_w\n    final_w=final_w/final_w.sum()\n\n    X=np.column_stack([data[f\"{a}__{t}\"].to_numpy(float) for a in ARMS])\n    p=X@final_w\n    dino=X[:,0]\n    pseudo_y=data[t].to_numpy(float)\n    gold_y=data[f\"gold__{t}\"].to_numpy(float)\n    non_gold=(~data[\"is_gold\"]).to_numpy()\n    gold_rows=data[\"is_gold\"].to_numpy()\n\n    pa=auc_soft(pseudo_y[non_gold],p[non_gold])\n    p0=auc_soft(pseudo_y[non_gold],dino[non_gold])\n    ga=auc_binary(gold_y[gold_rows],p[gold_rows])\n    g0=auc_binary(gold_y[gold_rows],dino[gold_rows])\n\n    gold_idx=np.flatnonzero(gold_rows & np.isfinite(gold_y))\n    boot_gains=[]\n    if len(gold_idx)>0:\n        yg=gold_y[gold_idx]\n        pos=gold_idx[yg==1]\n        neg=gold_idx[yg==0]\n        if len(pos)>0 and len(neg)>0:\n            for _ in range(GOLD_BOOTSTRAPS):\n                bi=np.concatenate([\n                    rng.choice(pos,size=len(pos),replace=True),\n                    rng.choice(neg,size=len(neg),replace=True)\n                ])\n                cg=auc_binary(gold_y[bi],p[bi])\n                bg=auc_binary(gold_y[bi],dino[bi])\n                if np.isfinite(cg) and np.isfinite(bg):\n                    boot_gains.append(cg-bg)\n\n    final_rows.append({\n        \"target\":t,\n        **{a:float(final_w[i]) for i,a in enumerate(ARMS)},\n        \"raw_best_dino\":float(raw_w[0]),\n        \"raw_best_a5\":float(raw_w[1]),\n        \"raw_best_rad\":float(raw_w[2]),\n        \"raw_best_b3\":float(raw_w[3]),\n        \"gold_auc_final\":ga,\n        \"gold_auc_dino\":g0,\n        \"gold_gain\":ga-g0 if np.isfinite(ga) and np.isfinite(g0) else np.nan,\n        \"pseudo_auc_final\":pa,\n        \"pseudo_auc_dino\":p0,\n        \"pseudo_gain\":pa-p0 if np.isfinite(pa) and np.isfinite(p0) else np.nan,\n        \"gold_bootstrap_gain_median\":float(np.median(boot_gains)) if boot_gains else np.nan,\n        \"gold_bootstrap_p_gain_gt0\":float(np.mean(np.asarray(boot_gains)>0)) if boot_gains else np.nan,\n    })\n\nfinal_df=pd.DataFrame(final_rows)\nfinal_df.to_csv(OUT/\"recommended_target_weights.csv\",index=False)\npd.concat(top_rows,ignore_index=True).to_csv(OUT/\"top10_candidates_per_target.csv\",index=False)\n\nweight_json={\n    \"seed\":SEED,\n    \"arms\":ARMS,\n    \"grid_step\":GRID_STEP,\n    \"constraints\":{\"dino_min\":DINO_MIN,\"a5_max\":A5_MAX,\"rad_max\":RAD_MAX,\"b3_max\":B3_MAX},\n    \"score_weights\":{\n        \"gold\":GOLD_SCORE_WEIGHT,\n        \"pseudo\":PSEUDO_SCORE_WEIGHT,\n        \"fold_nonloss\":FOLD_NONLOSS_WEIGHT,\n        \"fold_gain\":FOLD_GAIN_WEIGHT,\n        \"dino_prior_penalty\":DINO_PRIOR_PENALTY\n    },\n    \"target_shrink_to_global\":TARGET_SHRINK_TO_GLOBAL,\n    \"global_weights\":{a:float(global_w[i]) for i,a in enumerate(ARMS)},\n    \"target_weights\":{\n        r[\"target\"]:{a:float(r[a]) for a in ARMS}\n        for _,r in final_df.iterrows()\n    }\n}\n(OUT/\"target_weights.json\").write_text(json.dumps(weight_json,indent=2))\n\npy_lines=[\"TARGET_WEIGHTS = {\"]\nfor _,r in final_df.iterrows():\n    py_lines.append(\n        f'    {r[\"target\"]!r}: {{'\n        +\", \".join(f'{a!r}: {float(r[a]):.6f}' for a in ARMS)\n        +\"},\"\n    )\npy_lines.append(\"}\")\n(OUT/\"target_weights.py\").write_text(\"\\n\".join(py_lines)+\"\\n\")\n\nsummary=pd.DataFrame({\n    \"metric\":[\"macro Gold AUC\",\"macro pseudo AUC\"],\n    \"DINO\":[final_df[\"gold_auc_dino\"].mean(),final_df[\"pseudo_auc_dino\"].mean()],\n    \"OOF generalized blend\":[final_df[\"gold_auc_final\"].mean(),final_df[\"pseudo_auc_final\"].mean()]\n})\nsummary.to_csv(OUT/\"summary.csv\",index=False)\n\nprint(\"\\nRecommended target-specific weights after shrinkage:\")\nshow=final_df[[\"target\",\"dino\",\"a5\",\"rad\",\"b3\",\"gold_gain\",\"pseudo_gain\",\"gold_bootstrap_p_gain_gt0\"]].copy()\nprint(show.round(4).to_string(index=False))\n\nprint(\"\\nMacro comparison:\")\nprint(summary.round(4).to_string(index=False))\nprint(f\"\\nSaved -> {OUT}\")\nprint(\"Key files:\")\nprint(\"  target_weights.json\")\nprint(\"  target_weights.py\")\nprint(\"  recommended_target_weights.csv\")\nprint(\"  arm_rank_correlations.csv\")\nprint(\"  global_candidate_ranking.csv\")\nprint(\"  top10_candidates_per_target.csv\")\n","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","trusted":true,"execution":{"iopub.status.busy":"2026-08-16T02:24:31.332334Z","iopub.execute_input":"2026-08-16T02:24:31.332675Z","iopub.status.idle":"2026-08-16T02:26:10.387034Z","shell.execute_reply.started":"2026-08-16T02:24:31.332649Z","shell.execute_reply":"2026-08-16T02:26:10.385698Z"}},"outputs":[],"execution_count":null}]}