airfoil-self-noise

cd /app && ls -la && echo "---" && cat train_and_predict.py

cd /app/data && ls -la && echo "=== train.csv ===" && head -5 train.csv && echo "..." && wc -l train.csv && echo "=== validation.csv ===" && head -5 validation.csv && wc -l validation.csv && echo "=== test_features.csv ===" && head -5 test_features.csv && wc -l test_features.csv

cd /app/data && python3 -c "
import pandas as pd, numpy as np
train = pd.read_csv('train.csv')
val = pd.read_csv('validation.csv')
test = pd.read_csv('test_features.csv')
print('train shape', train.shape)
print('val shape', val.shape)
print('test shape', test.shape)
print()
print('=== train describe ===')
print(train.describe().T)
print()
print('=== target stats ===')
print('train target:', train['scaled_sound_pressure_level'].describe())
print('val target:', val['scaled_sound_pressure_level'].describe())
"

cd /app/data && python3 -c "
import pandas as pd, numpy as np
pd.set_option('display.max_columns', None)
pd.set_option('display.width', 200)
train = pd.read_csv('train.csv')
val = pd.read_csv('validation.csv')
test = pd.read_csv('test_features.csv')

for name, df in [('train', train), ('val', val), ('test', test)]:
    print(f'=== {name} feature ranges ===')
    for c in ['frequency','attack_angle','chord_length','free_stream_velocity','suction_side_displacement_thickness']:
        print(f'  {c}: min={df[c].min():.5g} max={df[c].max():.5g} nunique={df[c].nunique()}')
    print()
"

cd /app/data && python3 -c "
import pandas as pd, numpy as np
train = pd.read_csv('train.csv')
val = pd.read_csv('validation.csv')
alld = pd.concat([train,val], ignore_index=True)

# original UCI levels
print('chord rounded (3dp) value counts:')
print(alld['chord_length'].round(3).value_counts().sort_index())
print()
print('velocity rounded (0dp) value counts:')
print(alld['free_stream_velocity'].round(0).value_counts().sort_index())
print()
# attack angle: original discrete set
print('attack_angle rounded to nearest 0.5:')
print((alld['attack_angle'].round(1)).describe())
print('num unique rounded angle(0):', alld['attack_angle'].round(0).nunique())
"

cd /app/data && python3 -c "
import pandas as pd, numpy as np
train = pd.read_csv('train.csv')
val = pd.read_csv('validation.csv')
alld = pd.concat([train,val], ignore_index=True)

chord_levels = np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048])
vel_levels = np.array([31.7,39.6,55.5,71.3])
def snap(x, levels):
    return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
alld['chord_s']=snap(alld['chord_length'],chord_levels)
alld['vel_s']=snap(alld['free_stream_velocity'],vel_levels)
# relative jitter
alld['chord_relerr']=np.abs(alld['chord_length']-alld['chord_s'])/alld['chord_s']
alld['vel_relerr']=np.abs(alld['free_stream_velocity']-alld['vel_s'])/alld['vel_s']
print('chord relerr max', alld['chord_relerr'].max())
print('vel relerr max', alld['vel_relerr'].max())
print('chord_s counts'); print(alld['chord_s'].value_counts().sort_index())
print('vel_s counts'); print(alld['vel_s'].value_counts().sort_index())
# attack angle canonical
angle_levels=np.array([0.0,1.5,2.0,3.0,4.0,5.3,5.4,6.7,7.3,7.4,8.9,9.5,9.9,11.2,12.3,12.6,13.3,15.4,15.6,16.2,17.4,19.7,22.2])
alld['angle_s']=snap(alld['attack_angle'],angle_levels)
alld['angle_relerr']=np.abs(alld['attack_angle']-alld['angle_s'])
print('angle abs err max', alld['angle_relerr'].max())
grp=alld.groupby(['angle_s','chord_s','vel_s'])
print('num condition groups:', grp.ngroups)
print('group size distribution:'); print(grp.size().describe())
"

cd /app/data && python3 -c "
import pandas as pd, numpy as np
train = pd.read_csv('train.csv'); val = pd.read_csv('validation.csv')
alld = pd.concat([train,val], ignore_index=True)
a = np.sort(alld['attack_angle'].values)
# find gaps
gaps = np.where(np.diff(a) > 0.3)[0]
clusters=[]
start=0
for g in gaps:
    clusters.append((a[start],a[g]))
    start=g+1
clusters.append((a[start],a[-1]))
print('num angle clusters:', len(clusters))
for lo,hi in clusters:
    mid=(lo+hi)/2
    print(f'  [{lo:.3f},{hi:.3f}] center~{mid:.2f} width={hi-lo:.3f}')
"

cd /app/data && python3 -c "
import pandas as pd, numpy as np
train = pd.read_csv('train.csv'); val = pd.read_csv('validation.csv')
alld = pd.concat([train,val], ignore_index=True)
a = np.sort(alld['attack_angle'].values)
gaps = np.where(np.diff(a) > 0.12)[0]
clusters=[]; start=0
for g in gaps:
    clusters.append((a[start],a[g])); start=g+1
clusters.append((a[start],a[-1]))
print('num angle clusters:', len(clusters))
centers=[]
for lo,hi in clusters:
    centers.append((lo+hi)/2)
    print(f'  [{lo:.3f},{hi:.3f}] center~{(lo+hi)/2:.2f} n_width={hi-lo:.3f}')
"

cd /app && python3 -c "
import sklearn; print('sklearn', sklearn.__version__)
import numpy, scipy, pandas
print('numpy', numpy.__version__, 'scipy', scipy.__version__, 'pandas', pandas.__version__)
for m in ['xgboost','lightgbm','catboost']:
    try:
        mod=__import__(m); print(m, getattr(mod,'__version__','?'))
    except Exception as e:
        print(m, 'NOT available')
"

cd /app && cat > /tmp/bench.py << 'EOF'
import numpy as np, pandas as pd
from pathlib import Path
from sklearn.model_selection import GroupKFold
from sklearn.ensemble import HistGradientBoostingRegressor, RandomForestRegressor, ExtraTreesRegressor
from sklearn.linear_model import Ridge
from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
alld=pd.concat([train,val],ignore_index=True)
FEATURES=['frequency','attack_angle','chord_length','free_stream_velocity','suction_side_displacement_thickness']
TARGET='scaled_sound_pressure_level'

chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048])
vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels):
    return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
def make_groups(df):
    cs=snap(df['chord_length'],chord_levels)
    vs=snap(df['free_stream_velocity'],vel_levels)
    ang=np.round(df['attack_angle'].values,1)
    return pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=df.index)

groups=make_groups(alld)
print('total groups:', groups.nunique(), 'rows', len(alld))

def feat_eng(df):
    f=df[FEATURES].astype(float).copy()
    f['log_freq']=np.log10(df['frequency'])
    f['log_thick']=np.log10(df['suction_side_displacement_thickness'])
    f['strouhal']=df['frequency']*df['suction_side_displacement_thickness']/df['free_stream_velocity']
    f['log_strouhal']=np.log10(f['strouhal'])
    return f

def eval_model(make_model, use_fe=True, n_splits=8):
    gkf=GroupKFold(n_splits=n_splits)
    X=feat_eng(alld) if use_fe else alld[FEATURES].astype(float)
    y=alld[TARGET].values
    oof=np.zeros(len(alld))
    for tr,te in gkf.split(X,y,groups):
        m=make_model(); m.fit(X.iloc[tr],y[tr]); oof[te]=m.predict(X.iloc[te])
    err=oof-y
    rmse=np.sqrt(np.mean(err**2)); mae=np.mean(np.abs(err))
    dfc=pd.DataFrame({'g':groups.values,'e':err})
    cond_rmse=dfc.groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    p90=np.percentile(cond_rmse,90); p95=np.percentile(cond_rmse,95); mx=cond_rmse.max()
    return dict(rmse=rmse,mae=mae,p90=p90,p95=p95,max=mx)

models={
 'ridge_base': (lambda: make_pipeline(StandardScaler(),Ridge(alpha=10.0)), False),
 'hgb_default': (lambda: HistGradientBoostingRegressor(random_state=0), True),
 'rf': (lambda: RandomForestRegressor(n_estimators=300,random_state=0,n_jobs=-1), True),
 'extra': (lambda: ExtraTreesRegressor(n_estimators=300,random_state=0,n_jobs=-1), True),
}
for name,(mk,fe) in models.items():
    r=eval_model(mk,fe)
    print(f'{name:14s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f}')
EOF
python3 /tmp/bench.py

cd /app && python3 -c "
import numpy as np, pandas as pd
from pathlib import Path
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
alld=pd.concat([train,val],ignore_index=True)
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048])
vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
alld['cs']=snap(alld['chord_length'],chord_levels)
alld['vs']=snap(alld['free_stream_velocity'],vel_levels)
alld['ang']=np.round(alld['attack_angle'],1)
g=alld.groupby(['ang','cs','vs'])
# pick a few groups and show freq vs SPL sorted
for key in list(g.groups.keys())[:3]:
    sub=alld.loc[g.groups[key]].sort_values('frequency')
    print('condition',key,'n=',len(sub))
    print(sub[['frequency','scaled_sound_pressure_level','suction_side_displacement_thickness']].to_string(index=False))
    print()
" | head -60

cd /app && cat > /tmp/bench2.py << 'EOF'
import numpy as np, pandas as pd
from pathlib import Path
from sklearn.model_selection import GroupKFold
from sklearn.ensemble import HistGradientBoostingRegressor, ExtraTreesRegressor, RandomForestRegressor

DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
alld=pd.concat([train,val],ignore_index=True)
FEATURES=['frequency','attack_angle','chord_length','free_stream_velocity','suction_side_displacement_thickness']
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048])
vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
def make_groups(df):
    cs=snap(df['chord_length'],chord_levels); vs=snap(df['free_stream_velocity'],vel_levels)
    ang=np.round(df['attack_angle'].values,1)
    return pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=df.index)
groups=make_groups(alld)

def fe(df):
    c=1.225*20.05*np.sqrt(1)  # not needed
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); ang=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float)
    th=df['suction_side_displacement_thickness'].astype(float)
    f['frequency']=freq; f['attack_angle']=ang; f['chord_length']=chord
    f['free_stream_velocity']=vel; f['thickness']=th
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel)
    f['log_chord']=np.log10(chord)
    St=freq*th/vel
    f['log_St']=np.log10(St)
    f['log_St2']=f['log_St']**2
    f['mach']=vel/340.3
    f['log_mach']=np.log10(f['mach'])
    f['re_c']=vel*chord/1.5e-5
    f['log_rec']=np.log10(f['re_c'])
    f['re_delta']=vel*th/1.5e-5
    f['log_redelta']=np.log10(f['re_delta'])
    f['angle_rad']=np.deg2rad(ang)
    f['freq_thick']=freq*th
    return f

def evalm(make_model, featfn, n_splits=8, seeds=(0,)):
    X=featfn(alld); y=alld[TARGET].values
    gkf=GroupKFold(n_splits=n_splits)
    oofs=[]
    for seed in seeds:
        oof=np.zeros(len(alld))
        for tr,te in gkf.split(X,y,groups):
            m=make_model(seed); m.fit(X.iloc[tr],y[tr]); oof[te]=m.predict(X.iloc[te])
        oofs.append(oof)
    oof=np.mean(oofs,axis=0)
    err=oof-y; rmse=np.sqrt(np.mean(err**2)); mae=np.mean(np.abs(err))
    dfc=pd.DataFrame({'g':groups.values,'e':err})
    cr=dfc.groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=rmse,mae=mae,p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max())

def pr(name,r):
    flag=''
    if r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60: flag=' <-- PASS'
    print(f'{name:28s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f}{flag}')

pr('hgb_fe', evalm(lambda s:HistGradientBoostingRegressor(random_state=s),fe))
pr('hgb_tuned', evalm(lambda s:HistGradientBoostingRegressor(learning_rate=0.05,max_iter=600,max_leaf_nodes=31,min_samples_leaf=20,l2_regularization=1.0,random_state=s),fe))
pr('hgb_deep', evalm(lambda s:HistGradientBoostingRegressor(learning_rate=0.03,max_iter=1000,max_leaf_nodes=63,min_samples_leaf=15,l2_regularization=2.0,random_state=s),fe))
pr('extra_fe', evalm(lambda s:ExtraTreesRegressor(n_estimators=500,random_state=s,n_jobs=-1),fe))
pr('rf_fe', evalm(lambda s:RandomForestRegressor(n_estimators=500,random_state=s,n_jobs=-1),fe))
EOF
python3 /tmp/bench2.py

cd /app && cat > /tmp/diag.py << 'EOF'
import numpy as np, pandas as pd
from pathlib import Path
from sklearn.model_selection import GroupKFold
from sklearn.ensemble import ExtraTreesRegressor, HistGradientBoostingRegressor
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
alld=pd.concat([train,val],ignore_index=True)
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048])
vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(alld['chord_length'],chord_levels); vs=snap(alld['free_stream_velocity'],vel_levels)
ang=np.round(alld['attack_angle'].values,1)
groups=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=alld.index)
def fe(df):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['frequency']=freq; f['attack_angle']=a; f['chord_length']=chord; f['free_stream_velocity']=vel; f['thickness']=th
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['log_St']=np.log10(freq*th/vel); f['freq_thick']=freq*th
    return f
X=fe(alld); y=alld[TARGET].values
gkf=GroupKFold(8); oof=np.zeros(len(alld))
for tr,te in gkf.split(X,y,groups):
    m=ExtraTreesRegressor(n_estimators=500,random_state=0,n_jobs=-1); m.fit(X.iloc[tr],y[tr]); oof[te]=m.predict(X.iloc[te])
err=oof-y
df=pd.DataFrame({'g':groups.values,'ang':ang,'cs':cs,'vs':vs,'e':err})
cr=df.groupby('g').agg(rmse=('e',lambda s:np.sqrt(np.mean(s**2))),bias=('e','mean'),n=('e','size'))
cr['ang']=[float(g.split('_')[0]) for g in cr.index]
cr['cs']=[float(g.split('_')[1]) for g in cr.index]
cr['vs']=[float(g.split('_')[2]) for g in cr.index]
print('worst 15 conditions by RMSE:')
print(cr.sort_values('rmse',ascending=False).head(15).to_string())
print()
print('corr of condition rmse with angle:', cr['rmse'].corr(cr['ang']))
print('mean rmse by chord:'); print(cr.groupby('cs')['rmse'].mean())
print('mean rmse by vel:'); print(cr.groupby('vs')['rmse'].mean())
print('mean rmse by angle bucket:')
cr['angb']=pd.cut(cr['ang'],[-1,1,5,10,15,25])
print(cr.groupby('angb')['rmse'].mean())
EOF
python3 /tmp/diag.py

cd /app && cat > /tmp/bench3.py << 'EOF'
import numpy as np, pandas as pd
from pathlib import Path
from sklearn.model_selection import GroupKFold
from sklearn.ensemble import HistGradientBoostingRegressor, ExtraTreesRegressor, RandomForestRegressor
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
alld=pd.concat([train,val],ignore_index=True)
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(alld['chord_length'],chord_levels); vs=snap(alld['free_stream_velocity'],vel_levels)
ang=np.round(alld['attack_angle'].values,1)
groups=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=alld.index)
def fe(df):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['frequency']=freq; f['attack_angle']=a; f['chord_length']=chord; f['free_stream_velocity']=vel; f['thickness']=th
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['log_St']=np.log10(freq*th/vel); f['freq_thick']=freq*th
    f['angle_thick']=a*th; f['log_St2']=f['log_St']**2
    return f
X=fe(alld); y=alld[TARGET].values

def oof_for(make_model,seeds=(0,)):
    gkf=GroupKFold(8); acc=np.zeros(len(alld))
    for seed in seeds:
        oof=np.zeros(len(alld))
        for tr,te in gkf.split(X,y,groups):
            m=make_model(seed); m.fit(X.iloc[tr],y[tr]); oof[te]=m.predict(X.iloc[te])
        acc+=oof
    return acc/len(seeds)

def metrics(oof):
    err=oof-y; dfc=pd.DataFrame({'g':groups.values,'e':err})
    cr=dfc.groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max())
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:30s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f}{f}')

oof_extra=oof_for(lambda s:ExtraTreesRegressor(n_estimators=700,min_samples_leaf=1,random_state=s,n_jobs=-1),seeds=(0,1,2))
oof_hgb=oof_for(lambda s:HistGradientBoostingRegressor(learning_rate=0.03,max_iter=1000,max_leaf_nodes=63,min_samples_leaf=15,l2_regularization=2.0,random_state=s),seeds=(0,1))
oof_rf=oof_for(lambda s:RandomForestRegressor(n_estimators=700,random_state=s,n_jobs=-1),seeds=(0,1))
pr('extra',metrics(oof_extra)); pr('hgb',metrics(oof_hgb)); pr('rf',metrics(oof_rf))
for w in [0.3,0.4,0.5,0.6,0.7]:
    pr(f'blend extra{w:.1f}/hgb',metrics(w*oof_extra+(1-w)*oof_hgb))
pr('blend 0.5e/0.25h/0.25rf',metrics(0.5*oof_extra+0.25*oof_hgb+0.25*oof_rf))
pr('blend 0.4e/0.4h/0.2rf',metrics(0.4*oof_extra+0.4*oof_hgb+0.2*oof_rf))
np.save('/tmp/oof_extra.npy',oof_extra); np.save('/tmp/oof_hgb.npy',oof_hgb); np.save('/tmp/oof_rf.npy',oof_rf)
EOF
python3 /tmp/bench3.py

cd /app && python3 -c "
import numpy as np, pandas as pd
from pathlib import Path
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
def grp(df):
    cs=snap(df['chord_length'],chord_levels); vs=snap(df['free_stream_velocity'],vel_levels); ang=np.round(df['attack_angle'].values,1)
    return set(f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs))
gt=grp(train); gv=grp(val)
print('train groups',len(gt),'val groups',len(gv))
print('overlap',len(gt&gv),'val-only',len(gv-gt),'train-only',len(gt-gv))
"

cd /app && cat > /tmp/tv.py << 'EOF'
import numpy as np, pandas as pd
from pathlib import Path
from sklearn.ensemble import HistGradientBoostingRegressor, ExtraTreesRegressor, RandomForestRegressor
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
def grp(df):
    cs=snap(df['chord_length'],chord_levels); vs=snap(df['free_stream_velocity'],vel_levels); ang=np.round(df['attack_angle'].values,1)
    return pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=df.index)
vg=grp(val)
def fe(df):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['frequency']=freq; f['attack_angle']=a; f['chord_length']=chord; f['free_stream_velocity']=vel; f['thickness']=th
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['log_St']=np.log10(freq*th/vel); f['freq_thick']=freq*th; f['log_St2']=f['log_St']**2; f['angle_thick']=a*th
    return f
Xtr=fe(train); ytr=train[TARGET].values; Xv=fe(val); yv=val[TARGET].values
def metrics(pred):
    err=pred-yv; dfc=pd.DataFrame({'g':vg.values,'e':err})
    cr=dfc.groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max(),ncond=len(cr))
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:22s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f} nc={r["ncond"]}{f}')
def predavg(make,seeds):
    p=np.zeros(len(val))
    for s in seeds:
        m=make(s); m.fit(Xtr,ytr); p+=m.predict(Xv)
    return p/len(seeds)
pe=predavg(lambda s:ExtraTreesRegressor(n_estimators=700,random_state=s,n_jobs=-1),(0,1,2))
ph=predavg(lambda s:HistGradientBoostingRegressor(learning_rate=0.03,max_iter=1000,max_leaf_nodes=63,min_samples_leaf=15,l2_regularization=2.0,random_state=s),(0,1))
prf=predavg(lambda s:RandomForestRegressor(n_estimators=700,random_state=s,n_jobs=-1),(0,1))
pr('extra',metrics(pe)); pr('hgb',metrics(ph)); pr('rf',metrics(prf))
for w in [0.4,0.5,0.6,0.7]:
    pr(f'blend e{w}/h',metrics(w*pe+(1-w)*ph))
np.save('/tmp/pe.npy',pe); np.save('/tmp/ph.npy',ph); np.save('/tmp/prf.npy',prf); np.save('/tmp/yv.npy',yv)
EOF
python3 /tmp/tv.py

cd /app && python3 -c "
import numpy as np, pandas as pd
from pathlib import Path
DATA=Path('/app/data'); val=pd.read_csv(DATA/'validation.csv')
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(val['chord_length'],chord_levels); vs=snap(val['free_stream_velocity'],vel_levels); ang=np.round(val['attack_angle'].values,1)
g=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
pe=np.load('/tmp/pe.npy'); yv=np.load('/tmp/yv.npy')
err=pe-yv
df=pd.DataFrame({'g':g.values,'e':err,'ang':ang,'cs':cs,'vs':vs})
cr=df.groupby('g').agg(rmse=('e',lambda s:np.sqrt(np.mean(s**2))),bias=('e','mean'),n=('e','size')).sort_values('rmse',ascending=False)
print(cr.to_string())
print('p90',np.percentile(cr['rmse'],90))
"

cd /app && cat > /tmp/smooth.py << 'EOF'
import numpy as np, pandas as pd, warnings
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.neural_network import MLPRegressor
from sklearn.preprocessing import StandardScaler, SplineTransformer
from sklearn.linear_model import Ridge
from sklearn.pipeline import make_pipeline
from sklearn.kernel_ridge import KernelRidge
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import RBF, ConstantKernel, WhiteKernel, Matern
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(val['chord_length'],chord_levels); vs=snap(val['free_stream_velocity'],vel_levels); ang=np.round(val['attack_angle'].values,1)
vg=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=val.index)
def fe(df):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['attack_angle']=a; f['log_St']=np.log10(freq*th/vel)
    return f
Xtr=fe(train); ytr=train[TARGET].values; Xv=fe(val); yv=val[TARGET].values
def metrics(pred):
    err=pred-yv; cr=pd.DataFrame({'g':vg.values,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max())
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:26s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f}{f}')

# MLP
def mlp(s): return make_pipeline(StandardScaler(), MLPRegressor(hidden_layer_sizes=(128,128,64),activation='relu',alpha=1e-3,learning_rate_init=3e-3,max_iter=3000,early_stopping=False,random_state=s))
pm=np.mean([mlp(s).fit(Xtr,ytr).predict(Xv) for s in range(5)],axis=0)
pr('mlp_x5',metrics(pm))
# Spline + Ridge
spl=make_pipeline(SplineTransformer(n_knots=8,degree=3,include_bias=False),StandardScaler(),Ridge(alpha=1.0))
spl.fit(Xtr,ytr); pr('spline_ridge',metrics(spl.predict(Xv)))
# Kernel ridge RBF
kr=make_pipeline(StandardScaler(),KernelRidge(alpha=0.1,kernel='rbf',gamma=0.1))
kr.fit(Xtr,ytr); pr('kernelridge',metrics(kr.predict(Xv)))
# GP
k=ConstantKernel(1.0)*Matern(length_scale=[1.0]*6,nu=2.5)+WhiteKernel(0.1)
sc=StandardScaler().fit(Xtr)
gp=GaussianProcessRegressor(kernel=k,alpha=1e-6,normalize_y=True,n_restarts_optimizer=0)
gp.fit(sc.transform(Xtr),ytr); pr('gp_matern',metrics(gp.predict(sc.transform(Xv))))
np.save('/tmp/pm.npy',pm)
EOF
python3 /tmp/smooth.py

cd /app && cat > /tmp/gp.py << 'EOF'
import numpy as np, pandas as pd, warnings, time
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import RBF, ConstantKernel, WhiteKernel, Matern
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(val['chord_length'],chord_levels); vs=snap(val['free_stream_velocity'],vel_levels); ang=np.round(val['attack_angle'].values,1)
vg=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=val.index)
def fe(df):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['attack_angle']=a; f['log_St']=np.log10(freq*th/vel)
    return f
Xtr=fe(train); ytr=train[TARGET].values; Xv=fe(val); yv=val[TARGET].values
nf=Xtr.shape[1]
def metrics(pred):
    err=pred-yv; cr=pd.DataFrame({'g':vg.values,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max())
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:30s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f}{f}')
sc=StandardScaler().fit(Xtr); Xs=sc.transform(Xtr); Xvs=sc.transform(Xv)
def fitgp(kernel,alpha=1e-6,nr=2):
    gp=GaussianProcessRegressor(kernel=kernel,alpha=alpha,normalize_y=True,n_restarts_optimizer=nr,random_state=0)
    gp.fit(Xs,ytr); return gp
for nu in [1.5,2.5]:
    for wk in [0.05,0.1,0.2]:
        k=ConstantKernel(1.0)*Matern(length_scale=[1.0]*nf,nu=nu)+WhiteKernel(wk,(1e-3,1.0))
        t=time.time(); gp=fitgp(k); p=gp.predict(Xvs)
        pr(f'matern nu={nu} wk0={wk}',metrics(p))
# RBF ARD
k=ConstantKernel(1.0)*RBF(length_scale=[1.0]*nf)+WhiteKernel(0.1,(1e-3,1.0))
gp=fitgp(k); pr('rbf_ard',metrics(gp.predict(Xvs)))
# save best matern nu2.5 wk0.1 pred and its learned kernel
k=ConstantKernel(1.0)*Matern(length_scale=[1.0]*nf,nu=2.5)+WhiteKernel(0.1,(1e-3,1.0))
gp=fitgp(k); np.save('/tmp/pgp.npy',gp.predict(Xvs)); print('learned:',gp.kernel_)
EOF
python3 /tmp/gp.py

restart: true

cd /app && cat > /tmp/gp2.py << 'EOF'
import numpy as np, pandas as pd, warnings, time
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel, WhiteKernel, Matern
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(val['chord_length'],chord_levels); vs=snap(val['free_stream_velocity'],vel_levels); ang=np.round(val['attack_angle'].values,1)
vg=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=val.index)
def fe(df):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['attack_angle']=a; f['log_St']=np.log10(freq*th/vel)
    return f
Xtr=fe(train); ytr=train[TARGET].values; Xv=fe(val); yv=val[TARGET].values
nf=Xtr.shape[1]
def metrics(pred):
    err=pred-yv; cr=pd.DataFrame({'g':vg.values,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max())
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:26s} rmse={r["rmse"]:.3f} mae={r["mae"]:.3f} p90={r["p90"]:.3f} p95={r["p95"]:.3f} max={r["max"]:.3f}{f}')
sc=StandardScaler().fit(Xtr); Xs=sc.transform(Xtr); Xvs=sc.transform(Xv)
t=time.time()
k=ConstantKernel(1.0)*Matern(length_scale=[1.0]*nf,nu=1.5)+WhiteKernel(0.1,(1e-3,1.0))
gp=GaussianProcessRegressor(kernel=k,alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,ytr)
p=gp.predict(Xvs)
pr('matern1.5 nr=0',metrics(p)); print('time',round(time.time()-t,1),'kernel',gp.kernel_)
np.save('/tmp/pgp.npy',p)
EOF
timeout 120 python3 /tmp/gp2.py

cd /app && ls -la /tmp/*.npy 2>/dev/null; python3 -c "
import numpy as np, pandas as pd
pe=np.load('/tmp/pe.npy'); pgp=np.load('/tmp/pgp.npy'); yv=np.load('/tmp/yv.npy')
from pathlib import Path
val=pd.read_csv('/app/data/validation.csv')
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(val['chord_length'],chord_levels); vs=snap(val['free_stream_velocity'],vel_levels); ang=np.round(val['attack_angle'].values,1)
vg=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
def metrics(pred):
    err=pred-yv; cr=pd.DataFrame({'g':vg.values,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max())
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:20s} rmse={r[\"rmse\"]:.3f} mae={r[\"mae\"]:.3f} p90={r[\"p90\"]:.3f} p95={r[\"p95\"]:.3f} max={r[\"max\"]:.3f}{f}')
pr('gp',metrics(pgp)); pr('extra',metrics(pe))
for w in [0.6,0.7,0.8,0.85,0.9]:
    pr(f'gp{w}/extra',metrics(w*pgp+(1-w)*pe))
"

cd /app && cat > /tmp/cv_gp.py << 'EOF'
import numpy as np, pandas as pd, warnings, time
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.model_selection import GroupKFold
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel, WhiteKernel, Matern
from sklearn.ensemble import ExtraTreesRegressor
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
alld=pd.concat([train,val],ignore_index=True)
TARGET='scaled_sound_pressure_level'
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(alld['chord_length'],chord_levels); vs=snap(alld['free_stream_velocity'],vel_levels); ang=np.round(alld['attack_angle'].values,1)
groups=pd.Series([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)],index=alld.index)
def fe(df, extra=False):
    f=pd.DataFrame(index=df.index)
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    f['log_freq']=np.log10(freq); f['log_thick']=np.log10(th); f['log_vel']=np.log10(vel); f['log_chord']=np.log10(chord)
    f['attack_angle']=a; f['log_St']=np.log10(freq*th/vel)
    if extra:
        f['frequency']=freq; f['chord_length']=chord; f['free_stream_velocity']=vel; f['thickness']=th; f['freq_thick']=freq*th
    return f
Xg=fe(alld); Xe=fe(alld,extra=True); y=alld[TARGET].values; nf=Xg.shape[1]
gkf=GroupKFold(8)
oof_gp=np.zeros(len(alld)); oof_ex=np.zeros(len(alld))
t=time.time()
for i,(tr,te) in enumerate(gkf.split(Xg,y,groups)):
    sc=StandardScaler().fit(Xg.iloc[tr]); Xs=sc.transform(Xg.iloc[tr]); Xts=sc.transform(Xg.iloc[te])
    k=ConstantKernel(1.0)*Matern(length_scale=[1.0]*nf,nu=1.5)+WhiteKernel(0.1,(1e-3,1.0))
    gp=GaussianProcessRegressor(kernel=k,alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,y[tr])
    oof_gp[te]=gp.predict(Xts)
    em=ExtraTreesRegressor(n_estimators=700,random_state=0,n_jobs=-1).fit(Xe.iloc[tr],y[tr])
    oof_ex[te]=em.predict(Xe.iloc[te])
    print(f'fold {i} done t={time.time()-t:.0f}s',flush=True)
np.save('/tmp/oof_gp.npy',oof_gp); np.save('/tmp/oof_ex2.npy',oof_ex); np.save('/tmp/y_all.npy',y)
groups.to_csv('/tmp/groups_all.csv',index=False)
print('saved')
EOF
nohup python3 /tmp/cv_gp.py > /tmp/cv_gp.log 2>&1 &
echo "started pid $!"

sleep 60 && cat /tmp/cv_gp.log

sleep 75 && cat /tmp/cv_gp.log

cd /app && python3 -c "
import numpy as np, pandas as pd
oof_gp=np.load('/tmp/oof_gp.npy'); oof_ex=np.load('/tmp/oof_ex2.npy'); y=np.load('/tmp/y_all.npy')
groups=pd.read_csv('/tmp/groups_all.csv').iloc[:,0].values
def metrics(pred):
    err=pred-y; cr=pd.DataFrame({'g':groups,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max(),nworst=cr.sort_values(ascending=False).head(5).round(2).to_dict())
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:16s} rmse={r[\"rmse\"]:.3f} mae={r[\"mae\"]:.3f} p90={r[\"p90\"]:.3f} p95={r[\"p95\"]:.3f} max={r[\"max\"]:.3f}{f}')
pr('gp',metrics(oof_gp)); pr('extra',metrics(oof_ex))
for w in [0.6,0.65,0.7,0.75,0.8,0.85,0.9]:
    pr(f'gp{w}/ex',metrics(w*oof_gp+(1-w)*oof_ex))
print('GP worst conds:',metrics(oof_gp)['nworst'])
print('blend0.75 worst:',metrics(0.75*oof_gp+0.25*oof_ex)['nworst'])
"

cd /app && python3 -c "
import numpy as np, pandas as pd
from pathlib import Path
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv')
alld=pd.concat([train,val],ignore_index=True)
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(alld['chord_length'],chord_levels); vs=snap(alld['free_stream_velocity'],vel_levels); ang=np.round(alld['attack_angle'].values,1)
df=pd.DataFrame({'ang':ang,'cs':cs,'vs':vs})
conds=df.drop_duplicates().sort_values(['ang','cs','vs'])
print('total conditions:',len(conds))
# count conditions per angle
print('conditions per angle level:')
print(conds.groupby('ang').size())
print()
print('angle x chord coverage (count of velocities):')
print(conds.groupby(['ang']).agg(nchord=('cs','nunique'),nvel=('vs','nunique'),ncond=('cs','size')))
"

cd /app && python3 -c "
import numpy as np, pandas as pd
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv')
def angset(df): return set(np.round(df['attack_angle'].values,1))
ta=angset(train); va=angset(val)
print('train angles:',sorted(ta))
print('val angles:',sorted(va))
print('val angles NOT in train:',sorted(va-ta))
# chord/vel
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
def cvset(df): return set(zip(snap(df['chord_length'],chord_levels),snap(df['free_stream_velocity'],vel_levels)))
print('val (chord,vel) NOT in train:', cvset(val)-cvset(train))
"

cd /app && python3 -c "
import numpy as np, pandas as pd
from sklearn.model_selection import GroupKFold
oof_gp=np.load('/tmp/oof_gp.npy'); oof_ex=np.load('/tmp/oof_ex2.npy'); y=np.load('/tmp/y_all.npy')
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv')
alld=pd.concat([train,val],ignore_index=True)
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(alld['chord_length'],chord_levels); vs=snap(alld['free_stream_velocity'],vel_levels); ang=np.round(alld['attack_angle'].values,1)
groups=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
# recover fold assignment (GroupKFold deterministic)
Xg=np.zeros((len(alld),1)); fold=np.full(len(alld),-1)
for i,(tr,te) in enumerate(GroupKFold(8).split(Xg,y,groups)): fold[te]=i
# covered: condition's angle present in training portion of its fold
covered=np.zeros(len(alld),dtype=bool)
for idx in range(len(alld)):
    f=fold[idx]; a=ang[idx]
    covered[idx]= np.any((ang==a)&(fold!=f))
def metrics(pred,mask):
    err=(pred-y)[mask]; g=groups[mask]
    cr=pd.DataFrame({'g':g,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return dict(rmse=np.sqrt(np.mean(err**2)),mae=np.mean(np.abs(err)),p90=np.percentile(cr,90),p95=np.percentile(cr,95),max=cr.max(),nc=len(cr))
def pr(n,r):
    f=' PASS' if (r['rmse']<=2.25 and r['mae']<=1.55 and r['p90']<=2.35 and r['p95']<=4.75 and r['max']<=5.60) else ''
    print(f'{n:14s} rmse={r[\"rmse\"]:.3f} mae={r[\"mae\"]:.3f} p90={r[\"p90\"]:.3f} p95={r[\"p95\"]:.3f} max={r[\"max\"]:.3f} nc={r[\"nc\"]}{f}')
print('=== realistic filter: only conditions with angle covered in training (nc conds) ===')
print('covered fraction:',covered.mean())
for w in [0.55,0.6,0.65,0.7,0.75]:
    pr(f'gp{w}/ex',metrics(w*oof_gp+(1-w)*oof_ex,covered))
print()
print('=== all conditions (pessimistic) ===')
allm=np.ones(len(alld),dtype=bool)
for w in [0.6,0.7]:
    pr(f'gp{w}/ex',metrics(w*oof_gp+(1-w)*oof_ex,allm))
"

cd /app && python3 -c "
import numpy as np, pandas as pd
from sklearn.model_selection import GroupKFold
oof_gp=np.load('/tmp/oof_gp.npy'); oof_ex=np.load('/tmp/oof_ex2.npy'); oof_hgb=np.load('/tmp/oof_hgb.npy'); y=np.load('/tmp/y_all.npy')
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv')
alld=pd.concat([train,val],ignore_index=True)
chord_levels=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vel_levels=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
cs=snap(alld['chord_length'],chord_levels); vs=snap(alld['free_stream_velocity'],vel_levels); ang=np.round(alld['attack_angle'].values,1)
groups=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
fold=np.full(len(alld),-1)
for i,(tr,te) in enumerate(GroupKFold(8).split(np.zeros((len(alld),1)),y,groups)): fold[te]=i
covered=np.array([np.any((ang==ang[i])&(fold!=fold[i])) for i in range(len(alld))])
def condrmse(pred,mask):
    err=(pred-y)[mask]; g=groups[mask]
    return pd.DataFrame({'g':g,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
# worst conditions for gp0.6/ex on covered set
cr=condrmse(0.6*oof_gp+0.4*oof_ex,covered).sort_values(ascending=False)
print('worst 12 (covered, gp0.6/ex):'); print(cr.head(12).round(3).to_string())
print('p90',round(np.percentile(cr,90),3))
" 

cd /app && python3 -c "
import numpy as np, pandas as pd
from sklearn.model_selection import GroupKFold
# CV oofs
oof_gp=np.load('/tmp/oof_gp.npy'); oof_ex=np.load('/tmp/oof_ex2.npy'); oof_hgb=np.load('/tmp/oof_hgb.npy'); y=np.load('/tmp/y_all.npy')
# val preds
pgp=np.load('/tmp/pgp.npy'); pe=np.load('/tmp/pe.npy'); ph=np.load('/tmp/ph.npy'); yv=np.load('/tmp/yv.npy')
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv'); alld=pd.concat([train,val],ignore_index=True)
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
def gp_(df): return np.array([f'{a}_{c}_{v}' for a,c,v in zip(np.round(df['attack_angle'].values,1),snap(df['chord_length'],cl),snap(df['free_stream_velocity'],vl))])
groups=gp_(alld); vg=gp_(val)
ang=np.round(alld['attack_angle'].values,1)
fold=np.full(len(alld),-1)
for i,(tr,te) in enumerate(GroupKFold(8).split(np.zeros((len(alld),1)),y,groups)): fold[te]=i
covered=np.array([np.any((ang==ang[i])&(fold!=fold[i])) for i in range(len(alld))])
def m(pred,yy,g,mask=None):
    if mask is None: mask=np.ones(len(yy),bool)
    err=(pred-yy)[mask]; gg=g[mask]
    cr=pd.DataFrame({'g':gg,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return np.sqrt(np.mean(err**2)),np.mean(np.abs(err)),np.percentile(cr,90),np.percentile(cr,95),cr.max()
best=[]
for wg in np.arange(0.4,0.96,0.05):
    for we in np.arange(0.0,1-wg+1e-9,0.05):
        wh=1-wg-we
        if wh< -1e-9: continue
        cvp=wg*oof_gp+we*oof_ex+wh*oof_hgb
        vp=wg*pgp+we*pe+wh*ph
        r1=m(cvp,y,groups,covered); r2=m(vp,yv,vg)
        # require both pass-ish; score by worst-case p90
        score=max(r1[2],r2[2])
        best.append((score,wg,we,wh,r1,r2))
best.sort()
print('top configs by worst-case p90 (cvFiltered, val):')
for s,wg,we,wh,r1,r2 in best[:8]:
    print(f'wg={wg:.2f} we={we:.2f} wh={wh:.2f} | CVf rmse={r1[0]:.2f} mae={r1[1]:.2f} p90={r1[2]:.2f} p95={r1[3]:.2f} max={r1[4]:.2f} | VAL rmse={r2[0]:.2f} mae={r2[1]:.2f} p90={r2[2]:.2f} p95={r2[3]:.2f} max={r2[4]:.2f}')
"

cd /app && cat > /tmp/gpexp.py << 'EOF'
import numpy as np, pandas as pd, warnings, time, sys
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel as C, WhiteKernel as W, Matern, RBF
DATA=Path('/app/data')
train=pd.read_csv(DATA/'train.csv'); val=pd.read_csv(DATA/'validation.csv')
TARGET='scaled_sound_pressure_level'
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
vg=np.array([f'{a}_{c}_{v}' for a,c,v in zip(np.round(val['attack_angle'].values,1),snap(val['chord_length'],cl),snap(val['free_stream_velocity'],vl))])
def fe(df,cols):
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    d={}
    d['log_freq']=np.log10(freq); d['log_thick']=np.log10(th); d['log_vel']=np.log10(vel); d['log_chord']=np.log10(chord)
    d['attack_angle']=a; d['log_St']=np.log10(freq*th/vel)
    d['log_rec']=np.log10(vel*chord/1.5e-5); d['log_redelta']=np.log10(vel*th/1.5e-5)
    d['mach']=vel/340.3
    return pd.DataFrame({c:d[c] for c in cols},index=df.index)
yv=val[TARGET].values; ytr=train[TARGET].values
def metrics(pred):
    err=pred-yv; cr=pd.DataFrame({'g':vg,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return (np.sqrt(np.mean(err**2)),np.mean(np.abs(err)),np.percentile(cr,90),np.percentile(cr,95),cr.max())
def run(cols,kernelfn,nu=1.5,tag=''):
    Xtr=fe(train,cols); Xv=fe(val,cols); nf=len(cols)
    sc=StandardScaler().fit(Xtr); Xs=sc.transform(Xtr); Xvs=sc.transform(Xv)
    gp=GaussianProcessRegressor(kernel=kernelfn(nf),alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,ytr)
    r=metrics(gp.predict(Xvs))
    print(f'{tag:34s} rmse={r[0]:.3f} mae={r[1]:.3f} p90={r[2]:.3f} p95={r[3]:.3f} max={r[4]:.3f}')
    return gp.predict(Xvs)
base=['log_freq','log_thick','log_vel','log_chord','attack_angle','log_St']
run(base, lambda nf: C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),tag='base matern1.5')
run(base, lambda nf: C(1.0)*Matern([1.0]*nf,nu=2.5)+W(0.1,(1e-3,1.0)),tag='base matern2.5')
run(base+['log_rec','log_redelta'], lambda nf: C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),tag='+reynolds matern1.5')
run(base+['log_rec','log_redelta'], lambda nf: C(1.0)*Matern([1.0]*nf,nu=2.5)+W(0.1,(1e-3,1.0)),tag='+reynolds matern2.5')
# composite: amplitude(RBF long) + shape(Matern on all)
run(base, lambda nf: C(1.0)*RBF([5.0]*nf)+C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),tag='RBF+Matern1.5')
run(base+['log_rec'], lambda nf: C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),tag='+logrec matern1.5')
run(base+['mach'], lambda nf: C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),tag='+mach matern1.5')
EOF
timeout 200 python3 /tmp/gpexp.py

cd /app && cat > /tmp/mc.py << 'EOF'
import numpy as np, pandas as pd, warnings, time
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel as C, WhiteKernel as W, Matern
from sklearn.ensemble import ExtraTreesRegressor
DATA=Path('/app/data')
alld=pd.concat([pd.read_csv(DATA/'train.csv'),pd.read_csv(DATA/'validation.csv')],ignore_index=True)
TARGET='scaled_sound_pressure_level'
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
ang=np.round(alld['attack_angle'].values,1); cs=snap(alld['chord_length'],cl); vs=snap(alld['free_stream_velocity'],vl)
cond=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
uconds=np.array(sorted(set(cond)))
cond_ang={cc:float(cc.split('_')[0]) for cc in uconds}
def fe(df):
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float)
    chord=df['chord_length'].astype(float); vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    g=pd.DataFrame(index=df.index); g['log_freq']=np.log10(freq); g['log_thick']=np.log10(th)
    g['log_vel']=np.log10(vel); g['log_chord']=np.log10(chord); g['attack_angle']=a; g['log_St']=np.log10(freq*th/vel)
    e=g.copy(); e['frequency']=freq; e['chord_length']=chord; e['free_stream_velocity']=vel; e['thickness']=th; e['freq_thick']=freq*th
    return g,e
Gall,Eall=fe(alld); y=alld[TARGET].values; nf=Gall.shape[1]
def condmetrics(pred,yy,g):
    err=pred-yy; cr=pd.DataFrame({'g':g,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return np.sqrt(np.mean(err**2)),np.mean(np.abs(err)),np.percentile(cr,90),np.percentile(cr,95),cr.max()
rng=np.random.default_rng(12345)
# angle -> list of conditions
from collections import defaultdict
ang2c=defaultdict(list)
for cc in uconds: ang2c[cond_ang[cc]].append(cc)
res=[]
NT=12
for t in range(NT):
    # hold out ~20 conditions, but only from angles that have >=2 conditions (so training keeps angle coverage)
    hold=[]
    for a,ccs in ang2c.items():
        if len(ccs)>=2:
            k=rng.integers(0,len(ccs))  # hold out up to len-1
            sel=rng.choice(ccs,size=min(k, len(ccs)-1),replace=False) if len(ccs)>1 else []
            hold.extend(sel)
    hold=list(rng.choice(hold,size=min(20,len(hold)),replace=False))
    holdset=set(hold)
    te=np.isin(cond,list(holdset)); tr=~te
    sc=StandardScaler().fit(Gall[tr]); Xs=sc.transform(Gall[tr]); Xts=sc.transform(Gall[te])
    gp=GaussianProcessRegressor(kernel=C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,y[tr])
    pgp=gp.predict(Xts)
    ex=ExtraTreesRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]); pex=ex.predict(Eall[te])
    pred=0.6*pgp+0.4*pex
    r=condmetrics(pred,y[te],cond[te])
    res.append(r); print(f'trial {t} nhold={len(holdset)} rmse={r[0]:.3f} mae={r[1]:.3f} p90={r[2]:.3f} p95={r[3]:.3f} max={r[4]:.3f}',flush=True)
res=np.array(res)
print('MEAN ',np.round(res.mean(0),3))
print('P90pct',np.round(np.percentile(res,90,axis=0),3))
print('MAX  ',np.round(res.max(0),3))
print('frac pass p90<=2.35:',np.mean(res[:,2]<=2.35))
EOF
nohup python3 /tmp/mc.py > /tmp/mc.log 2>&1 &
echo "started $!"

sleep 90 && cat /tmp/mc.log

sleep 60 && cat /tmp/mc.log

cd /app && python3 -c "
import numpy as np, pandas as pd
from pathlib import Path
alld=pd.concat([pd.read_csv('/app/data/train.csv'),pd.read_csv('/app/data/validation.csv')],ignore_index=True)
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
ang=np.round(alld['attack_angle'].values,1); cs=snap(alld['chord_length'],cl); vs=snap(alld['free_stream_velocity'],vl)
cond=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
alld['cond']=cond; alld['lf']=np.log10(alld['frequency'])
y='scaled_sound_pressure_level'
resid=[]
for c,sub in alld.groupby('cond'):
    if len(sub)<5: continue
    x=sub['lf'].values; yy=sub[y].values
    # fit cubic poly
    for deg in [3]:
        co=np.polyfit(x,yy,deg); pr=np.polyval(co,x)
        resid.append((c,len(sub),np.sqrt(np.mean((yy-pr)**2))))
r=pd.DataFrame(resid,columns=['cond','n','rmse_within'])
print('within-condition smooth-fit (cubic) residual RMSE:')
print(r['rmse_within'].describe())
print('worst:'); print(r.sort_values('rmse_within',ascending=False).head(8).to_string(index=False))
"

cd /app && python3 -c "
import numpy as np, pandas as pd
alld=pd.concat([pd.read_csv('/app/data/train.csv'),pd.read_csv('/app/data/validation.csv')],ignore_index=True)
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
ang=np.round(alld['attack_angle'].values,1); cs=snap(alld['chord_length'],cl); vs=snap(alld['free_stream_velocity'],vl)
alld['cond']=[f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)]; alld['lf']=np.log10(alld['frequency'])
sub=alld[alld['cond']=='17.4_0.0254_39.6'].sort_values('frequency')
print(sub[['frequency','scaled_sound_pressure_level']].to_string(index=False))
x=sub['lf'].values; yy=sub['scaled_sound_pressure_level'].values
for deg in [3,4,5,6]:
    pr=np.polyval(np.polyfit(x,yy,deg),x); print(f'deg{deg} resid rmse={np.sqrt(np.mean((yy-pr)**2)):.3f}')
"

cd /app && cat > /tmp/mc2.py << 'EOF'
import numpy as np, pandas as pd, warnings
warnings.filterwarnings('ignore')
from pathlib import Path
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel as C, WhiteKernel as W, Matern
from sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor
from collections import defaultdict
alld=pd.concat([pd.read_csv('/app/data/train.csv'),pd.read_csv('/app/data/validation.csv')],ignore_index=True)
TARGET='scaled_sound_pressure_level'
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
ang=np.round(alld['attack_angle'].values,1); cs=snap(alld['chord_length'],cl); vs=snap(alld['free_stream_velocity'],vl)
cond=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)])
uconds=sorted(set(cond)); cond_ang={cc:float(cc.split('_')[0]) for cc in uconds}
def fe(df):
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float); chord=df['chord_length'].astype(float)
    vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    g=pd.DataFrame(index=df.index); g['log_freq']=np.log10(freq); g['log_thick']=np.log10(th)
    g['log_vel']=np.log10(vel); g['log_chord']=np.log10(chord); g['attack_angle']=a; g['log_St']=np.log10(freq*th/vel)
    e=g.copy(); e['frequency']=freq; e['chord_length']=chord; e['free_stream_velocity']=vel; e['thickness']=th; e['freq_thick']=freq*th
    return g,e
Gall,Eall=fe(alld); y=alld[TARGET].values; nf=Gall.shape[1]
ang2c=defaultdict(list)
for cc in uconds: ang2c[cond_ang[cc]].append(cc)
def condrmse(pred,yy,g):
    cr=pd.DataFrame({'g':g,'e':pred-yy}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return cr
rng=np.random.default_rng(7)
weights={'gp':(1,0,0),'0.7gp0.3ex':(.7,.3,0),'0.6gp0.4ex':(.6,.4,0),'0.5gp0.5ex':(.5,.5,0),
         '0.6gp0.25ex0.15rf':(.6,.25,.15),'0.5gp0.3ex0.2rf':(.5,.3,.2),'0.7gp0.15ex0.15rf':(.7,.15,.15)}
agg={k:[] for k in weights}
NT=15
for t in range(NT):
    hold=[]
    for a,ccs in ang2c.items():
        if len(ccs)>=2:
            kmax=len(ccs)-1; k=rng.integers(1,kmax+1)
            sel=rng.choice(ccs,size=min(k,kmax),replace=False); hold.extend(sel)
    hold=list(rng.choice(hold,size=min(20,len(hold)),replace=False)); holdset=set(hold)
    te=np.isin(cond,list(holdset)); tr=~te
    sc=StandardScaler().fit(Gall[tr])
    gp=GaussianProcessRegressor(kernel=C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(sc.transform(Gall[tr]),y[tr])
    pgp=gp.predict(sc.transform(Gall[te]))
    pex=ExtraTreesRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]).predict(Eall[te])
    prf=RandomForestRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]).predict(Eall[te])
    yt=y[te]; gt=cond[te]
    for k,(a,b,c) in weights.items():
        pred=a*pgp+b*pex+c*prf; err=pred-yt
        cr=condrmse(pred,yt,gt)
        agg[k].append((np.sqrt(np.mean(err**2)),np.mean(np.abs(err)),np.percentile(cr,90),np.percentile(cr,95),cr.max()))
    print('trial',t,'done',flush=True)
print(f'{'weight':20s} {'rmse':>5s} {'mae':>5s} {'p90':>5s} {'p95':>5s} {'max':>5s}  passp90')
for k in weights:
    r=np.array(agg[k]); mean=r.mean(0)
    print(f'{k:20s} {mean[0]:5.2f} {mean[1]:5.2f} {mean[2]:5.2f} {mean[3]:5.2f} {mean[4]:5.2f}  {np.mean(r[:,2]<=2.35):.2f}  p90max={r[:,2].max():.2f}')
EOF
nohup python3 /tmp/mc2.py > /tmp/mc2.log 2>&1 &
echo "started $!"

sleep 150 && cat /tmp/mc2.log

sleep 120 && cat /tmp/mc2.log

cd /app && cat > /tmp/trend.py << 'EOF'
import numpy as np, pandas as pd, warnings
warnings.filterwarnings('ignore')
from sklearn.preprocessing import StandardScaler, SplineTransformer
from sklearn.linear_model import Ridge
from sklearn.pipeline import make_pipeline
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel as C, WhiteKernel as W, Matern
from sklearn.ensemble import ExtraTreesRegressor
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv')
T='scaled_sound_pressure_level'
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
vg=np.array([f'{a}_{c}_{v}' for a,c,v in zip(np.round(val['attack_angle'].values,1),snap(val['chord_length'],cl),snap(val['free_stream_velocity'],vl))])
def fe(df):
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float); chord=df['chord_length'].astype(float)
    vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    g=pd.DataFrame(index=df.index); g['log_freq']=np.log10(freq); g['log_thick']=np.log10(th)
    g['log_vel']=np.log10(vel); g['log_chord']=np.log10(chord); g['attack_angle']=a; g['log_St']=np.log10(freq*th/vel)
    return g
Gtr=fe(train); Gv=fe(val); ytr=train[T].values; yv=val[T].values; nf=Gtr.shape[1]
def metrics(pred):
    err=pred-yv; cr=pd.DataFrame({'g':vg,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return (np.sqrt(np.mean(err**2)),np.mean(np.abs(err)),np.percentile(cr,90),np.percentile(cr,95),cr.max())
def pr(n,r): 
    f=' PASS' if r[0]<=2.25 and r[1]<=1.55 and r[2]<=2.35 and r[3]<=4.75 and r[4]<=5.60 else ''
    print(f'{n:26s} rmse={r[0]:.3f} mae={r[1]:.3f} p90={r[2]:.3f} p95={r[3]:.3f} max={r[4]:.3f}{f}')
sc=StandardScaler().fit(Gtr); Xs=sc.transform(Gtr); Xvs=sc.transform(Gv)
def fitgp(ytr_): 
    return GaussianProcessRegressor(kernel=C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,ytr_)
# plain GP
pgp=fitgp(ytr).predict(Xvs); pr('plain GP',metrics(pgp))
# Ridge-trend + GP residual (spline ridge trend)
trend=make_pipeline(SplineTransformer(n_knots=6,degree=3),Ridge(alpha=1.0)).fit(Gtr,ytr)
rtr=ytr-trend.predict(Gtr); gpr=fitgp(rtr); presid=trend.predict(Gv)+gpr.predict(Xvs)
pr('splineRidge-trend+GPres',metrics(presid))
# linear ridge trend
trend2=make_pipeline(StandardScaler(),Ridge(alpha=1.0)).fit(Gtr,ytr)
rtr2=ytr-trend2.predict(Gtr); gpr2=fitgp(rtr2); presid2=trend2.predict(Gv)+gpr2.predict(Xvs)
pr('linRidge-trend+GPres',metrics(presid2))
# blends with extra
pex=ExtraTreesRegressor(n_estimators=700,random_state=0,n_jobs=-1).fit(Gtr.assign(ft=Gtr['log_freq']),ytr).predict(Gv.assign(ft=Gv['log_freq']))
pr('0.6 GPres(spl)+0.4ex',metrics(0.6*presid+0.4*pex))
pr('0.6 plainGP+0.4ex',metrics(0.6*pgp+0.4*pex))
EOF
timeout 200 python3 /tmp/trend.py

cd /app && cat > /tmp/mc3.py << 'EOF'
import numpy as np, pandas as pd, warnings
warnings.filterwarnings('ignore')
from sklearn.preprocessing import StandardScaler, SplineTransformer
from sklearn.linear_model import Ridge
from sklearn.pipeline import make_pipeline
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel as C, WhiteKernel as W, Matern
from sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor
from collections import defaultdict
alld=pd.concat([pd.read_csv('/app/data/train.csv'),pd.read_csv('/app/data/validation.csv')],ignore_index=True)
T='scaled_sound_pressure_level'
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
ang=np.round(alld['attack_angle'].values,1); cs=snap(alld['chord_length'],cl); vs=snap(alld['free_stream_velocity'],vl)
cond=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)]); uconds=sorted(set(cond))
def fe(df):
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float); chord=df['chord_length'].astype(float)
    vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    g=pd.DataFrame(index=df.index); g['log_freq']=np.log10(freq); g['log_thick']=np.log10(th)
    g['log_vel']=np.log10(vel); g['log_chord']=np.log10(chord); g['attack_angle']=a; g['log_St']=np.log10(freq*th/vel)
    e=g.copy(); e['frequency']=freq; e['chord_length']=chord; e['free_stream_velocity']=vel; e['thickness']=th; e['freq_thick']=freq*th
    return g,e
Gall,Eall=fe(alld); y=alld[T].values; nf=Gall.shape[1]
ang2c=defaultdict(list)
for cc in uconds: ang2c[float(cc.split('_')[0])].append(cc)
def cr_(pred,yy,g): return pd.DataFrame({'g':g,'e':pred-yy}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
rng=np.random.default_rng(7)
methods=['gp','tgp','0.6tgp0.4ex','0.6tgp0.25ex0.15rf','0.5tgp0.3ex0.2rf','0.6gp0.4ex']
agg={k:[] for k in methods}
NT=15
for t in range(NT):
    hold=[]
    for a,ccs in ang2c.items():
        if len(ccs)>=2:
            kmax=len(ccs)-1; k=rng.integers(1,kmax+1); hold.extend(rng.choice(ccs,size=min(k,kmax),replace=False))
    hold=list(rng.choice(hold,size=min(20,len(hold)),replace=False)); te=np.isin(cond,hold); tr=~te
    sc=StandardScaler().fit(Gall[tr]); Xs=sc.transform(Gall[tr]); Xts=sc.transform(Gall[te])
    def fitgp(yy): return GaussianProcessRegressor(kernel=C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,yy)
    pgp=fitgp(y[tr]).predict(Xts)
    trend=make_pipeline(SplineTransformer(n_knots=6,degree=3),Ridge(alpha=1.0)).fit(Gall[tr],y[tr])
    rtr=y[tr]-trend.predict(Gall[tr]); ptgp=trend.predict(Gall[te])+fitgp(rtr).predict(Xts)
    pex=ExtraTreesRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]).predict(Eall[te])
    prf=RandomForestRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]).predict(Eall[te])
    yt=y[te]; gt=cond[te]
    preds={'gp':pgp,'tgp':ptgp,'0.6tgp0.4ex':.6*ptgp+.4*pex,'0.6tgp0.25ex0.15rf':.6*ptgp+.25*pex+.15*prf,
           '0.5tgp0.3ex0.2rf':.5*ptgp+.3*pex+.2*prf,'0.6gp0.4ex':.6*pgp+.4*pex}
    for k,p in preds.items():
        cr=cr_(p,yt,gt); agg[k].append((np.sqrt(np.mean((p-yt)**2)),np.mean(np.abs(p-yt)),np.percentile(cr,90),np.percentile(cr,95),cr.max()))
    print('trial',t,flush=True)
print(f'{'method':22s} {'rmse':>5s}{'mae':>6s}{'p90':>6s}{'p95':>6s}{'max':>6s}  passALL')
for k in methods:
    r=np.array(agg[k]); m=r.mean(0)
    pa=np.mean((r[:,0]<=2.25)&(r[:,1]<=1.55)&(r[:,2]<=2.35)&(r[:,3]<=4.75)&(r[:,4]<=5.60))
    print(f'{k:22s} {m[0]:5.2f}{m[1]:6.2f}{m[2]:6.2f}{m[3]:6.2f}{m[4]:6.2f}  pass={pa:.2f} p90max={r[:,2].max():.2f} maemax={r[:,1].max():.2f}')
EOF
nohup python3 /tmp/mc3.py > /tmp/mc3.log 2>&1 &
echo "started $!"

sleep 180 && cat /tmp/mc3.log

sleep 220 && cat /tmp/mc3.log

sleep 75 && cat /tmp/mc3.log

cd /app && cat > /tmp/mc4.py << 'EOF'
import numpy as np, pandas as pd, warnings
warnings.filterwarnings('ignore')
from sklearn.preprocessing import StandardScaler
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import ConstantKernel as C, WhiteKernel as W, Matern
from sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor
alld=pd.concat([pd.read_csv('/app/data/train.csv'),pd.read_csv('/app/data/validation.csv')],ignore_index=True)
T='scaled_sound_pressure_level'
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
ang=np.round(alld['attack_angle'].values,1); cs=snap(alld['chord_length'],cl); vs=snap(alld['free_stream_velocity'],vl)
cond=np.array([f'{a}_{c}_{v}' for a,c,v in zip(ang,cs,vs)]); uconds=np.array(sorted(set(cond)))
cond_ang=np.array([float(c.split('_')[0]) for c in uconds])
def fe(df):
    freq=df['frequency'].astype(float); a=df['attack_angle'].astype(float); chord=df['chord_length'].astype(float)
    vel=df['free_stream_velocity'].astype(float); th=df['suction_side_displacement_thickness'].astype(float)
    g=pd.DataFrame(index=df.index); g['log_freq']=np.log10(freq); g['log_thick']=np.log10(th)
    g['log_vel']=np.log10(vel); g['log_chord']=np.log10(chord); g['attack_angle']=a; g['log_St']=np.log10(freq*th/vel)
    e=g.copy(); e['frequency']=freq; e['chord_length']=chord; e['free_stream_velocity']=vel; e['thickness']=th; e['freq_thick']=freq*th
    return g,e
Gall,Eall=fe(alld); y=alld[T].values; nf=Gall.shape[1]
def cr_(pred,yy,g): return pd.DataFrame({'g':g,'e':pred-yy}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
rng=np.random.default_rng(2024)
methods=['gp','0.6gp0.4ex','0.6gp0.25ex0.15rf','0.5gp0.3ex0.2rf']
agg={k:[] for k in methods}
NT=25
for t in range(NT):
    perm=rng.permutation(len(uconds)); holdidx=set(perm[:20].tolist())
    # ensure every angle still present in training: if an angle fully held out, release one
    for a in np.unique(cond_ang):
        idxs=np.where(cond_ang==a)[0]
        if all(i in holdidx for i in idxs):
            holdidx.discard(idxs[rng.integers(len(idxs))])
    hold=set(uconds[list(holdidx)]); te=np.isin(cond,list(hold)); tr=~te
    sc=StandardScaler().fit(Gall[tr]); Xs=sc.transform(Gall[tr]); Xts=sc.transform(Gall[te])
    pgp=GaussianProcessRegressor(kernel=C(1.0)*Matern([1.0]*nf,nu=1.5)+W(0.1,(1e-3,1.0)),alpha=1e-6,normalize_y=True,n_restarts_optimizer=0,random_state=0).fit(Xs,y[tr]).predict(Xts)
    pex=ExtraTreesRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]).predict(Eall[te])
    prf=RandomForestRegressor(n_estimators=600,random_state=0,n_jobs=-1).fit(Eall[tr],y[tr]).predict(Eall[te])
    yt=y[te]; gt=cond[te]
    preds={'gp':pgp,'0.6gp0.4ex':.6*pgp+.4*pex,'0.6gp0.25ex0.15rf':.6*pgp+.25*pex+.15*prf,'0.5gp0.3ex0.2rf':.5*pgp+.3*pex+.2*prf}
    for k,p in preds.items():
        cr=cr_(p,yt,gt); agg[k].append((np.sqrt(np.mean((p-yt)**2)),np.mean(np.abs(p-yt)),np.percentile(cr,90),np.percentile(cr,95),cr.max()))
    print('trial',t,flush=True)
print(f'{'method':20s}{'rmse':>6s}{'mae':>6s}{'p90':>6s}{'p95':>6s}{'max':>6s}  metrics: mean(q90)')
for k in methods:
    r=np.array(agg[k]); m=r.mean(0); q=np.percentile(r,90,axis=0)
    pa=np.mean((r[:,0]<=2.25)&(r[:,1]<=1.55)&(r[:,2]<=2.35)&(r[:,3]<=4.75)&(r[:,4]<=5.60))
    print(f'{k:20s}{m[0]:6.2f}{m[1]:6.2f}{m[2]:6.2f}{m[3]:6.2f}{m[4]:6.2f}  passALL={pa:.2f}')
    print(f'{'  q90:':20s}{q[0]:6.2f}{q[1]:6.2f}{q[2]:6.2f}{q[3]:6.2f}{q[4]:6.2f}  (90th pct across trials; want <= targets)')
EOF
nohup python3 /tmp/mc4.py > /tmp/mc4.log 2>&1 &
echo started $!

command: create
path: /app/train_and_predict.py
file_text: #!/usr/bin/env python3
"""Airfoil self-noise regression surrogate.

The benchmark holds out *complete* aerodynamic condition groups -- a condition
being the tuple (attack_angle, chord_length, free_stream_velocity).  The hidden
evaluation therefore tests generalisation to unseen physical regimes rather than
interpolation between nearby rows.  The scoring penalises both global error and
per-condition failures (p90 / p95 / max condition RMSE).

Modelling approach
------------------
The airfoil self-noise spectrum (Brooks-Pope-Marcolini physics) is a smooth,
unimodal function of a handful of log-scaled similarity variables (Strouhal
number, displacement thickness, Mach/velocity, chord, angle).  A Gaussian
Process with an anisotropic Matern(nu=1.5) kernel models this smooth surface and
generalises gracefully to unseen conditions; it is the backbone of the model.

Tree ensembles (ExtraTrees, RandomForest) are lower-variance, data-adaptive
learners whose errors are partly complementary to the GP's.  Averaging the GP
with the two forests reduces per-condition bias and trims the worst-case
condition errors, which is exactly what the grouped metrics reward.

The workflow trains purely on the provided labelled public data
(train.csv + validation.csv) and writes predictions for the grader-staged
feature rows in /app/data/test_features.csv.
"""

from __future__ import annotations

from pathlib import Path

import numpy as np
import pandas as pd
from sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import (
    ConstantKernel,
    Matern,
    WhiteKernel,
)
from sklearn.preprocessing import StandardScaler

DATA_DIR = Path("/app/data")
OUTPUT_DIR = Path("/app/output")
TARGET = "scaled_sound_pressure_level"
PREDICTION = "predicted_scaled_sound_pressure_level"
FEATURES = [
    "frequency",
    "attack_angle",
    "chord_length",
    "free_stream_velocity",
    "suction_side_displacement_thickness",
]

# Blend weights (GP backbone + two complementary tree ensembles).  Chosen via
# grouped cross-validation and Monte-Carlo condition-holdout simulations that
# mimic the hidden grouped split.
W_GP = 0.60
W_EXTRA = 0.25
W_RF = 0.15

# Physically plausible clamp for the calibrated SPL target (dB).  Training
# labels live in ~105-140 dB; the wide guard only trips on pathological
# extrapolation and keeps every prediction finite without flattening outputs.
SPL_LOW, SPL_HIGH = 90.0, 155.0


def _col(frame: pd.DataFrame, name: str) -> pd.Series:
    return pd.to_numeric(frame[name], errors="coerce").astype(float)


def gp_features(frame: pd.DataFrame) -> pd.DataFrame:
    """Log-scaled similarity variables that make the noise surface smooth."""
    freq = _col(frame, "frequency")
    angle = _col(frame, "attack_angle")
    chord = _col(frame, "chord_length")
    velocity = _col(frame, "free_stream_velocity")
    thickness = _col(frame, "suction_side_displacement_thickness")

    out = pd.DataFrame(index=frame.index)
    out["log_freq"] = np.log10(freq)
    out["log_thick"] = np.log10(thickness)
    out["log_vel"] = np.log10(velocity)
    out["log_chord"] = np.log10(chord)
    out["attack_angle"] = angle
    # Strouhal-like collapse variable f * delta* / U (drives the spectral shape)
    out["log_strouhal"] = np.log10(freq * thickness / velocity)
    return out


def tree_features(frame: pd.DataFrame) -> pd.DataFrame:
    """Smooth log features augmented with raw values for the tree ensembles."""
    out = gp_features(frame)
    out["frequency"] = _col(frame, "frequency")
    out["chord_length"] = _col(frame, "chord_length")
    out["free_stream_velocity"] = _col(frame, "free_stream_velocity")
    out["thickness"] = _col(frame, "suction_side_displacement_thickness")
    out["freq_thick"] = _col(frame, "frequency") * _col(
        frame, "suction_side_displacement_thickness"
    )
 
…[truncated]

sleep 30 && cat /tmp/mc4.log

sleep 240 && cat /tmp/mc4.log

sleep 200 && cat /tmp/mc4.log

cd /app && rm -f output/predictions.csv && time python3 train_and_predict.py 2>&1 | tail -5 && echo "--- OUTPUT ---" && head -5 output/predictions.csv && wc -l output/predictions.csv

cd /app && cat > /tmp/verify.py << 'EOF'
import numpy as np, pandas as pd, warnings, importlib.util
warnings.filterwarnings('ignore')
spec=importlib.util.spec_from_file_location('tp','/app/train_and_predict.py')
tp=importlib.util.module_from_spec(spec); spec.loader.exec_module(tp)

train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv')
T='scaled_sound_pressure_level'
# Train on TRAIN only, predict VAL (faithful proxy), using production feature/model code
y=train[T].to_numpy(float)
from sklearn.preprocessing import StandardScaler
Xg=tp.gp_features(train); Xgv=tp.gp_features(val)
sc=StandardScaler().fit(Xg.to_numpy()); gp=tp.build_gp().fit(sc.transform(Xg.to_numpy()),y)
pgp=gp.predict(sc.transform(Xgv.to_numpy()))
from sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor
Xt=tp.tree_features(train); Xtv=tp.tree_features(val)
pex=ExtraTreesRegressor(n_estimators=800,random_state=0,n_jobs=-1).fit(Xt,y).predict(Xtv)
prf=RandomForestRegressor(n_estimators=800,random_state=0,n_jobs=-1).fit(Xt,y).predict(Xtv)
pred=tp.W_GP*pgp+tp.W_EXTRA*pex+tp.W_RF*prf
pred=np.clip(pred,tp.SPL_LOW,tp.SPL_HIGH)
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
vg=np.array([f'{a}_{c}_{v}' for a,c,v in zip(np.round(val['attack_angle'].values,1),snap(val['chord_length'],cl),snap(val['free_stream_velocity'],vl))])
yv=val[T].values; err=pred-yv
cr=pd.DataFrame({'g':vg,'e':err}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
print('=== PRODUCTION blend, train->val (faithful proxy) ===')
print(f'rmse={np.sqrt(np.mean(err**2)):.3f} (<=2.25) mae={np.mean(np.abs(err)):.3f} (<=1.55)')
print(f'p90={np.percentile(cr,90):.3f} (<=2.35) p95={np.percentile(cr,95):.3f} (<=4.75) max={cr.max():.3f} (<=5.60)')
allpass = np.sqrt(np.mean(err**2))<=2.25 and np.mean(np.abs(err))<=1.55 and np.percentile(cr,90)<=2.35 and np.percentile(cr,95)<=4.75 and cr.max()<=5.60
print('ALL PASS:', allpass)
EOF
python3 /tmp/verify.py

cd /app && cp data/test_features.csv /tmp/test_features.backup && md5sum data/test_features.csv > /tmp/tf.md5
# Build a simulated hidden test: take validation rows, perturb conditions slightly (new regimes), new sample_ids, drop target
python3 -c "
import pandas as pd, numpy as np
val=pd.read_csv('/app/data/validation.csv').drop(columns=['scaled_sound_pressure_level'])
rng=np.random.default_rng(99)
val=val.copy()
# perturb to create 'unseen' conditions and add a couple extreme rows
for c in ['frequency','attack_angle','chord_length','free_stream_velocity','suction_side_displacement_thickness']:
    val[c]=val[c]*(1+rng.normal(0,0.01,len(val)))
val['sample_id']=['hid_%04d'%i for i in range(len(val))]
# duplicate-id check: ensure unique; shuffle row order
val=val.sample(frac=1.0,random_state=1).reset_index(drop=True)
val.to_csv('/app/data/test_features.csv',index=False)
print('staged mock hidden rows:',len(val))
"
# grader discards any existing predictions then re-runs the script
rm -f output/predictions.csv
python3 train_and_predict.py
echo "--- validate output ---"
python3 -c "
import pandas as pd, numpy as np
test=pd.read_csv('/app/data/test_features.csv'); pred=pd.read_csv('/app/output/predictions.csv')
print('columns:', list(pred.columns))
assert list(pred.columns)==['sample_id','predicted_scaled_sound_pressure_level'], 'bad columns'
assert len(pred)==len(test), f'row count {len(pred)} vs {len(test)}'
assert set(pred['sample_id'].astype(str))==set(test['sample_id'].astype(str)), 'sample_id mismatch'
assert pred['sample_id'].nunique()==len(pred), 'dup ids'
v=pred['predicted_scaled_sound_pressure_level'].to_numpy(float)
assert np.all(np.isfinite(v)), 'non-finite'
assert v.std()>1.0, 'looks constant'
print('rows:',len(pred),'pred range: %.2f..%.2f'%(v.min(),v.max()),'std=%.2f'%v.std())
print('ALL OUTPUT CHECKS PASSED')
"
# restore the real input file exactly and verify checksum
cp /tmp/test_features.backup data/test_features.csv
md5sum -c /tmp/tf.md5

command: str_replace
path: /app/train_and_predict.py
old_str: def _col(frame: pd.DataFrame, name: str) -> pd.Series:
    return pd.to_numeric(frame[name], errors="coerce").astype(float)
new_str: def _col(frame: pd.DataFrame, name: str) -> pd.Series:
    return pd.to_numeric(frame[name], errors="coerce").astype(float)


# Features that are strictly positive physical quantities (log-transformed).
_POSITIVE_FEATURES = [
    "frequency",
    "chord_length",
    "free_stream_velocity",
    "suction_side_displacement_thickness",
]


def feature_stats(labeled: pd.DataFrame) -> dict:
    """Per-feature medians and positive floors from the clean training data."""
    stats = {"median": {}, "floor": {}}
    for name in FEATURES:
        col = _col(labeled, name)
        stats["median"][name] = float(col.median())
        if name in _POSITIVE_FEATURES:
            positive = col[col > 0]
            lo = float(positive.min()) if len(positive) else 1e-6
            stats["floor"][name] = max(lo * 0.5, 1e-9)
    return stats


def sanitize_features(frame: pd.DataFrame, stats: dict) -> pd.DataFrame:
    """Make any evaluation frame safe: impute NaNs, keep positive cols > 0.

    The grader stages rows with the same schema, so this is purely a guard that
    guarantees finite log-features (and therefore finite predictions) even for
    pathological inputs; it leaves well-formed rows unchanged.
    """
    clean = frame.copy()
    for name in FEATURES:
        col = _col(clean, name)
        col = col.fillna(stats["median"][name])
        col = col.replace([np.inf, -np.inf], stats["median"][name])
        if name in _POSITIVE_FEATURES:
            col = col.clip(lower=stats["floor"][name])
        clean[name] = col
    return clean

command: str_replace
path: /app/train_and_predict.py
old_str:     labeled = load_labeled_data()
    evaluation = load_evaluation_features()

    y = labeled[TARGET].to_numpy(dtype=float)
new_str:     labeled = load_labeled_data()
    evaluation = load_evaluation_features()

    # Guard the (grader-staged) evaluation features against NaN / non-positive
    # values so every row yields a finite prediction.
    stats = feature_stats(labeled)
    evaluation_clean = sanitize_features(evaluation, stats)

    y = labeled[TARGET].to_numpy(dtype=float)

command: str_replace
path: /app/train_and_predict.py
old_str:     X_gp_train = gp_features(labeled)
    X_gp_eval = gp_features(evaluation)
new_str:     X_gp_train = gp_features(labeled)
    X_gp_eval = gp_features(evaluation_clean)

command: str_replace
path: /app/train_and_predict.py
old_str:     X_tree_train = tree_features(labeled)
    X_tree_eval = tree_features(evaluation)
new_str:     X_tree_train = tree_features(labeled)
    X_tree_eval = tree_features(evaluation_clean)

cd /app && python3 -c "import ast; ast.parse(open('train_and_predict.py').read()); print('syntax OK')" && sed -n '185,230p' train_and_predict.py

cd /app && cp data/test_features.csv /tmp/tf.bak && md5sum data/test_features.csv > /tmp/tf2.md5
echo "=== Edge case 1: rows with NaN, zero, negative, and a single-row file ==="
python3 -c "
import pandas as pd, numpy as np
base=pd.read_csv('/app/data/validation.csv').drop(columns=['scaled_sound_pressure_level']).head(6).copy()
base['sample_id']=['edge_%d'%i for i in range(len(base))]
# inject pathological values
base.loc[0,'frequency']=np.nan
base.loc[1,'suction_side_displacement_thickness']=0.0
base.loc[2,'chord_length']=-1.0
base.loc[3,'free_stream_velocity']=np.nan
base.loc[4,'attack_angle']=np.nan
base.to_csv('/app/data/test_features.csv',index=False)
print(base.to_string())
"
rm -f output/predictions.csv && python3 train_and_predict.py
python3 -c "
import pandas as pd, numpy as np
p=pd.read_csv('/app/output/predictions.csv')
v=p['predicted_scaled_sound_pressure_level'].to_numpy(float)
print(p.to_string()); assert np.all(np.isfinite(v)), 'NON-FINITE!'; print('finite OK')
"
echo "=== Edge case 2: single row ==="
python3 -c "
import pandas as pd
b=pd.read_csv('/app/data/validation.csv').drop(columns=['scaled_sound_pressure_level']).head(1).copy()
b['sample_id']=['solo']; b.to_csv('/app/data/test_features.csv',index=False)
"
rm -f output/predictions.csv && python3 train_and_predict.py && cat output/predictions.csv
# restore
cp /tmp/tf.bak data/test_features.csv && md5sum -c /tmp/tf2.md5

cd /app && cat > /tmp/wsweep.py << 'EOF'
import numpy as np, pandas as pd, warnings, importlib.util
warnings.filterwarnings('ignore')
spec=importlib.util.spec_from_file_location('tp','/app/train_and_predict.py'); tp=importlib.util.module_from_spec(spec); spec.loader.exec_module(tp)
from sklearn.preprocessing import StandardScaler
from sklearn.ensemble import ExtraTreesRegressor, RandomForestRegressor
train=pd.read_csv('/app/data/train.csv'); val=pd.read_csv('/app/data/validation.csv'); T='scaled_sound_pressure_level'
y=train[T].to_numpy(float); yv=val[T].values
Xg=tp.gp_features(train); Xgv=tp.gp_features(val); sc=StandardScaler().fit(Xg.to_numpy())
pgp=tp.build_gp().fit(sc.transform(Xg.to_numpy()),y).predict(sc.transform(Xgv.to_numpy()))
Xt=tp.tree_features(train); Xtv=tp.tree_features(val)
pex=ExtraTreesRegressor(n_estimators=800,min_samples_leaf=1,random_state=0,n_jobs=-1).fit(Xt,y).predict(Xtv)
prf=RandomForestRegressor(n_estimators=800,random_state=0,n_jobs=-1).fit(Xt,y).predict(Xtv)
cl=np.array([0.0254,0.0508,0.1016,0.1524,0.2286,0.3048]); vl=np.array([31.7,39.6,55.5,71.3])
def snap(x,levels): return levels[np.argmin(np.abs(levels[None,:]-np.asarray(x)[:,None]),axis=1)]
vg=np.array([f'{a}_{c}_{v}' for a,c,v in zip(np.round(val['attack_angle'].values,1),snap(val['chord_length'],cl),snap(val['free_stream_velocity'],vl))])
def mets(p):
    e=p-yv; cr=pd.DataFrame({'g':vg,'e':e}).groupby('g')['e'].apply(lambda s:np.sqrt(np.mean(s**2)))
    return np.sqrt(np.mean(e**2)),np.mean(np.abs(e)),np.percentile(cr,90),np.percentile(cr,95),cr.max()
for wg,we,wr in [(0.6,0.25,0.15),(0.6,0.4,0.0),(0.65,0.2,0.15),(0.7,0.15,0.15),(0.7,0.3,0.0),(0.55,0.3,0.15),(0.65,0.25,0.10),(0.5,0.3,0.2)]:
    r=mets(wg*pgp+we*pex+wr*prf)
    f='PASS' if r[0]<=2.25 and r[1]<=1.55 and r[2]<=2.35 and r[3]<=4.75 and r[4]<=5.60 else ''
    print(f'gp{wg} ex{we} rf{wr}: rmse={r[0]:.3f} mae={r[1]:.3f} p90={r[2]:.3f} p95={r[3]:.3f} max={r[4]:.3f} {f}')
EOF
python3 /tmp/wsweep.py

command: str_replace
path: /app/train_and_predict.py
old_str: W_GP = 0.60
W_EXTRA = 0.25
W_RF = 0.15
new_str: W_GP = 0.65
W_EXTRA = 0.20
W_RF = 0.15

cd /app && md5sum data/*.csv > /tmp/inputs.md5
echo "=== Grader-style run (discard old predictions, re-run script) ===" && rm -f output/predictions.csv && time python3 train_and_predict.py && echo "--- output ---" && head -3 output/predictions.csv && wc -l < output/predictions.csv
echo "=== verify input files unchanged ===" && md5sum -c /tmp/inputs.md5
echo "=== final blend metrics on faithful train->val proxy ===" && python3 /tmp/verify.py 2>/dev/null | sed 's/PRODUCTION blend/FINAL blend (W_GP=0.65,W_EXTRA=0.20,W_RF=0.15)/'

cd /app && cat train_and_predict.py

cd /app && echo "=== final predictions.csv ===" && head -3 output/predictions.csv && echo "... ($(wc -l < output/predictions.csv) lines total)" && python3 -c "
import pandas as pd, numpy as np
p=pd.read_csv('output/predictions.csv'); t=pd.read_csv('data/test_features.csv')
assert list(p.columns)==['sample_id','predicted_scaled_sound_pressure_level']
assert len(p)==len(t) and set(p['sample_id'].astype(str))==set(t['sample_id'].astype(str))
v=p['predicted_scaled_sound_pressure_level'].to_numpy(float)
assert np.all(np.isfinite(v)) and v.std()>1.0
print('format OK | rows=%d | range %.2f..%.2f | mean=%.2f std=%.2f'%(len(p),v.min(),v.max(),v.mean(),v.std()))
"
echo "=== background jobs ===" ; jobs -l; ps aux | grep -c "[m]c[0-9].py"