{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.7.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":38257,"databundleVersionId":4319132,"sourceType":"competition"}],"dockerImageVersionId":30407,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# Model Error Analysis Example\n\n**The notebook illustrates drawing histograms of angular and azimuth errors using the weighted Line-fit method as an example.**\n\n## Unit system\n\nWe work in a system of units where coordinates are measured in kilometers and times in nanoseconds are multiplied by `0.299792458e-3`. As a result, the speed of light in vacuum is $c=1$. The conversion to these systems of units is done in the `get_sensors` and `prepare_batch` functions.","metadata":{}},{"cell_type":"markdown","source":"## Some Useful Functions","metadata":{}},{"cell_type":"code","source":"import os, gc, sys, time, datetime, math, random,  psutil\nimport numpy as np,  matplotlib.pyplot as plt, pandas as pd\nfrom   pathlib   import Path        \nfrom   tqdm.auto import tqdm\nimport pyarrow, pyarrow.parquet as pq     # read by chanks\n\nPATH = Path(\"/kaggle/input/icecube-neutrinos-in-deep-ice\")  # path to dataset\nfiles_trn = [item for item in (PATH  / \"train\").glob('*')]  # all train files\nprint(f\"{len(files_trn):3d} train files\")\n\n#===============================================================================\ndef info(text, pref=\"\", end=\"\\n\"):\n    \"\"\" Information about the progress of calculations (time and memory) \"\"\"\n    gc.collect()\n    ram, t = psutil.virtual_memory().used / 1024**3,  time.time()    \n    print(f\"{pref}{(t-info.beg)/60:5.1f}m[{t-info.last:+5.1f}s] {ram:6.3f}Gb > {text}\",end=end)\n    info.last = time.time(); \ninfo.beg = info.last = time.time()\n\n#-------------------------------------------------------------------------------\n\ndef get_sensors():\n    \"\"\" Get sensor positions \"\"\"    \n    df = pd.read_csv(PATH / \"sensor_geometry.csv\")      \n    df['line_id'] = df.sensor_id // 60 + 1                 # string id\n    df['core']    = (df.line_id > 78).astype(np.float32)   # sensor from DeepCore\n    df.x = ( df.x * 1e-3 ).astype(np.float32)              # distances in kilometers\n    df.y = ( df.y * 1e-3 ).astype(np.float32)\n    df.z = ( df.z * 1e-3 ).astype(np.float32)\n    \n    return df[['sensor_id', 'line_id', 'core', 'x', 'y', 'z']]\n\n\n#-------------------------------------------------------------------------------\n\ndef get_target_angles(batches=set(), verbose = True):\n    \"\"\" Get target angles for batches from  batches id from set batches \"\"\"    \n    res = None\n    file = pq.ParquetFile(PATH / \"train_meta.parquet\")    \n    for b in file.iter_batches(batch_size=200_000, columns=['event_id','batch_id','azimuth','zenith']):    \n        if len(batches) == 0:\n            break\n        batch_df = b.to_pandas()\n        batch_id = batch_df.batch_id[0]\n        if batch_id in batches:      \n            batch_df.event_id= batch_df.event_id.astype(np.int64)      \n            batch_df.azimuth = batch_df.azimuth.astype(np.float32)\n            batch_df.zenith  = batch_df.zenith.astype(np.float32)    \n            batch_df = batch_df[ ['event_id','batch_id','azimuth','zenith'] ]\n            res =  batch_df if res is None else res.append(batch_df)            \n            batches.remove(batch_id)            \n            if verbose: info(f\"load target for batch {batch_id};  {len(batches)} candies left\")\n    return res\n#-------------------------------------------------------------------------------\n\ndef prepare_batch(df, drop_aux = True, doms_agg = False, verbose=True):\n    \"\"\" \n    Preparing a loaded batch, shifting and normalizing times \n        * drop_aux - whether it is necessary to throw out auxiliary pulses (auxiliary==True)\n        * doms_agg - whether it is necessary to aggregate pulses on sensors\n        * verbose  - debug output\n    \"\"\"    \n    df['event_id'] = df.index.astype(np.int64)\n    df = df.reset_index(drop=True)  # sensor_id, t, charge, aux, event_id    \n    df.rename(columns={\"time\": \"t\", \"auxiliary\": \"aux\", 'charge': 'q'}, inplace=True)\n    df.q = df.q.astype(np.float32)\n\n    if drop_aux:\n        df = df[ ~df.aux ]\n    \n    if doms_agg:\n        df = df.groupby(['event_id', 'sensor_id']).agg(\n            aux = ( 'aux', \"mean\"),\n            q   = ( 'q',   \"sum\"),\n            t   = ( 't',   \"min\"),)\n        df = df.reset_index()    \n        \n    if verbose: info(f\"prepare_batch: loaded  {df.shape}\")\n        \n    times = df.groupby('event_id').agg( t_min = ('t', 'min') )\n    df = df.merge(times, left_on='event_id', right_index=True, how='left')\n    df.t = (( df.t - df.t_min ) * 0.299792458e-3 ).astype(np.float32)             \n    \n    if verbose: info(\"prepare_batch: shift_times\")    \n\n    return df[['event_id', 'sensor_id', 'aux', 'q', 't' ]]\n\n#-------------------------------------------------------------------------------\n\ndef agg_features(df):\n    \"\"\" Calculate some features aggregated for a given event \"\"\"\n    df = df.groupby('event_id').agg(\n        lines = ( 'line_id',   'nunique'),\n        doms  = ( 'sensor_id', 'nunique'),\n        core  = ( 'core',      'mean'),\n        pulses= ( 't',         'count'),\n        t     = ( 't',         'mean'),\n        q     = ( 'q',         'mean'),\n    )\n    return df\n#-------------------------------------------------------------------------------\n\ndef angles2vector(df):\n    \"\"\" Add unit vector components from (azimuth,zenith) to the DataFrame df \"\"\"\n    df['nx'] = np.sin(df.zenith) * np.cos(df.azimuth)\n    df['ny'] = np.sin(df.zenith) * np.sin(df.azimuth)\n    df['nz'] = np.cos(df.zenith) \n    return df\n\n#-------------------------------------------------------------------------------\n\ndef vector2angles(n, eps=1e-8):\n    \"\"\"  Get spherical angles of vector n: (B,3) \"\"\"                \n    n = n / (np.linalg.norm(n, axis=1, keepdims=True) + eps)    \n                                \n    azimuth = np.arctan2( n[:,1],  n[:,0])    \n    azimuth[azimuth < 0] += 2*np.pi\n                                \n    zenith = np.arccos( n[:,2].clip(-1,1) )                                \n    \n    return azimuth, zenith\n\n#-------------------------------------------------------------------------------\n\ndef angle_errors(n1, n2, eps=1e-8):\n    \"\"\" Calculate angles between two vectors:: n1,n2: (B,3) return: (B,) \"\"\"\n    n1 = n1 / (np.linalg.norm(n1, axis=1, keepdims=True) + eps)\n    n2 = n2 / (np.linalg.norm(n2, axis=1, keepdims=True) + eps)\n    \n    cos = (n1*n2).sum(axis=1)                     # angles between vectors\n    angle_err = np.arccos( cos.clip(-1,1) )    \n        \n    r1   =  n1[:,0]*n1[:,0] + n1[:,1]*n1[:,1]    # angles between vectors in (x,y)    \n    r2   =  n2[:,0]*n2[:,0] + n2[:,1]*n2[:,1]\n    cosX = (n1[:,0]*n2[:,0] + n1[:,1]*n2[:,1]) / (np.sqrt(r1*r2) + eps)    \n    azimuth_err = np.arccos( cosX.clip(-1,1) )\n                                \n    zerros = r1 < eps                            # azimuth angle not defined\n    azimuth_err[zerros] = np.random.random((len(n1[zerros]),))*np.pi\n    \n    zenith1  = np.arccos( n1[:,2].clip(-1,1) )\n    zenith2  = np.arccos( n2[:,2].clip(-1,1) )\n    zenith_err = zenith2 - zenith1    \n        \n    return angle_err, azimuth_err, zenith_err\n\n#-------------------------------------------------------------------------------\n\ndef hist(x, x_sub, pref=\"\", digits=1, w=12, h=3, bins=50, bins_sub=100):\n    \"\"\" Distribution of values of x and its subset x_sub (in a narrower range) \"\"\"\n    r = lambda x: '{x:.{digits}f}'.format(x=round(x,digits), digits=digits)\n    plt.figure(figsize=(w,h), facecolor ='w')         \n    plt.suptitle(f\"{pref}median={r(x.median())}; mean={r(x.mean())} ± {r(x.std())}  [min,max]=[{r(x.min())}, {r(x.max())}]; cnt={len(x)} ({100*len(x_sub)/len(x):.0f}%)\", fontsize=14)\n    plt.subplot(1,2,1)\n    plt.hist(x, bins=bins, log=True, color=\"lightblue\", ec=\"black\");  plt.grid(ls=\":\",alpha=1); plt.ylabel(\"log10(N)\")\n    plt.subplot(1,2,2)\n    plt.hist(x_sub, bins=bins_sub, density=True, color=\"lightblue\", ec=\"black\");    plt.grid(ls=\":\",alpha=1); plt.ylabel(\"Density\")\n    plt.show()","metadata":{"_uuid":"8f2839f25d086af736a60e9eeb907d3b93b6e0e5","_cell_guid":"b1076dfc-b9ad-4769-8c92-a6c4dae69d19","execution":{"iopub.status.busy":"2024-05-12T12:29:06.658934Z","iopub.execute_input":"2024-05-12T12:29:06.659687Z","iopub.status.idle":"2024-05-12T12:29:06.824746Z","shell.execute_reply.started":"2024-05-12T12:29:06.659626Z","shell.execute_reply":"2024-05-12T12:29:06.823225Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Load train batches and target angles","metadata":{}},{"cell_type":"code","source":"class CFG:\n    drop_aux    = True       # discard pulses with auxiliary == True; the method Line-fit needs\n    doms_agg    = False      # aggregate pulses on sensors\n    show_stat   = True       # draw histograms of values aggregated by events\n    batches     = range(1,6) # list of batch id-s    \n    \n    def get(end=\", \"):\n        return \"\".join([f\"{k}:{v}{end}\" for k,v in CFG.__dict__.items() if not k.startswith(\"__\") and k != \"get\"])\n\n#======================================================================================    \ninfo.beg = info.last = time.time()\ninfo(\"Begin\")\n\n# 1. load sensors\nsensors_df = get_sensors()   \nprint(\"sensor_df: \")\ndisplay(sensors_df.head(2))\n\n# 2. load pulses from batches\npulses_df = None              \nfor batch_id in tqdm(CFG.batches):\n    df = pd.read_parquet(PATH  / \"train\" / f\"batch_{batch_id}.parquet\")     \n    df = prepare_batch(df, drop_aux = CFG.drop_aux, doms_agg = CFG.doms_agg)\n    df = df.merge(sensors_df, left_on='sensor_id', right_on='sensor_id', how=\"left\")\n    pulses_df = df if pulses_df is None else pulses_df.append(df)\ndel df; \ndisplay(pulses_df.head(2))\n\n# 3. load target angles and get unit direction vector\ntrue_df = get_target_angles(batches=set( CFG.batches ))\ntrue_df = angles2vector(true_df)\ndisplay(true_df.head(2))\ninfo(\"loaded target\")\n        \ninfo(\"the End\")","metadata":{"execution":{"iopub.status.busy":"2024-05-12T12:29:06.833402Z","iopub.execute_input":"2024-05-12T12:29:06.833928Z","iopub.status.idle":"2024-05-12T12:30:14.957418Z","shell.execute_reply.started":"2024-05-12T12:29:06.833879Z","shell.execute_reply":"2024-05-12T12:30:14.955962Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Line-fit method\n\nThe velocity of the charge $\\mathbf{u}$ and the vector $\\mathbf{q}$ through which it passes at time $t=0$ are found by the following formulas:\n\n$$\n\\mathbf{u} = \\frac{\\langle\\mathbf{r}\\,t\\rangle - \\langle\\mathbf{r}\\rangle\\,\\langle t\\rangle}{\\langle t^2\\rangle - \\langle t\\rangle^2},~~~~~~~~~ \\mathbf{q} = \\langle\\mathbf{r}\\rangle - \\mathbf{u}\\,\\langle t\\rangle.\n$$\n\nThe averaging of $\\mathbf{r}=(x,y,z)$ coordinates and sensor response times can be done with weights:\n\n$$\n\\langle f\\rangle = \\sum_i w_i\\,f_i,~~~~~~~~~~~~~~\\sum_i w_i = 1.\n$$","metadata":{}},{"cell_type":"code","source":"def line_fit(df, lm = 5.5, eps = 1e-8):\n    \"\"\" Weighted Line-fit method \"\"\"\n    df['w']   = np.exp(-lm * df.t)      \n    \n    df['xw']  = df.x * df.w;     df['yw']  = df.y*df.w;       df['zw']  = df.z*df.w;  \n    df['xtw'] = df.x*df.t*df.w;  df['ytw'] = df.y*df.t*df.w;  df['ztw'] = df.z*df.t*df.w; \n    df['ttw'] = df.t*df.t*df.w;  df['tw']  = df.t*df.w; \n    agg = df.groupby([\"event_id\"]).agg(\n        xw  = ('xw', 'sum'),  yw  = ('yw', 'sum'),  zw  = ('zw', 'sum'),  tw = ('tw', 'sum'),\n        xtw = ('xtw','sum'),  ytw = ('ytw','sum'),  ztw = ('ztw','sum'),  ttw = ('ttw','sum'),\n        w   = ('w',  'sum')\n    )      \n    agg.xw  /= agg.w;  agg.yw  /= agg.w; agg.zw  /= agg.w;  agg.tw  /= agg.w;\n    agg.xtw /= agg.w;  agg.ytw /= agg.w; agg.ztw /= agg.w;  agg.ttw /= agg.w;\n    \n    agg['Dtw'] = agg.ttw - agg.tw*agg.tw\n    agg['ux'] = ( agg.xtw - agg.xw*agg.tw ) / ( agg.Dtw + eps);  agg['qx'] = agg.xw - agg.ux * agg.tw\n    agg['uy'] = ( agg.ytw - agg.yw*agg.tw ) / ( agg.Dtw + eps);  agg['qy'] = agg.yw - agg.uy * agg.tw\n    agg['uz'] = ( agg.ztw - agg.zw*agg.tw ) / ( agg.Dtw + eps);  agg['qz'] = agg.zw - agg.uz * agg.tw\n    return agg[['ux', 'uy', 'uz', 'qx', 'qy', 'qz' ]]\n    \n#--------------------------------------------------------------------------------------------\n\npred_df = line_fit(pulses_df, lm=5.5)\nu = pred_df[['ux','uy','uz']].to_numpy()\nn = true_df[['nx','ny','nz']].to_numpy()    \nang_err, az_err, ze_err = angle_errors(-u,n)                     # inverse vector!\ninfo(f\"mean={np.mean(ang_err):.3f}   medan={np.median(ang_err):.3f}   num={len(ang_err)}   min,max=[{np.min(ang_err):.3f},{np.max(ang_err):.3f}] {CFG.get()}\")\n","metadata":{"execution":{"iopub.status.busy":"2024-05-12T12:30:14.959384Z","iopub.execute_input":"2024-05-12T12:30:14.95984Z","iopub.status.idle":"2024-05-12T12:30:33.608037Z","shell.execute_reply.started":"2024-05-12T12:30:14.95979Z","shell.execute_reply":"2024-05-12T12:30:33.606159Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Error Analysis","metadata":{}},{"cell_type":"code","source":"def plot_hist_error(err, ax, bins):    \n    \"\"\" Plot a histogram of errors \"\"\"\n    #ax.set_facecolor(\"ivory\"); \n    ax.autoscale(tight=True)\n    p,_,_ = ax.hist(err, bins=bins, range=(0,np.pi), fc=\"lightblue\", density=True, alpha=0.8)   \n    ax.grid(ls=\":\"); ax.set_ylim([0, p.max()])\n    x = np.linspace(0, np.pi,bins) + (np.pi/bins)/2\n    ax.plot(x, np.cumsum(p)*(np.pi/bins), \":b\", lw=0.8 )\n    return p, x, np.pi/bins\n\n#-------------------------------------------------------------------------------\n\ndef plot_angle_error(err, ax, bins = 100):    \n    \"\"\" Plot a histogram of angle errors and calculate the statistics and the share w of 'bad examples' \"\"\"\n    p, x, dx = plot_hist_error(err, ax, bins)\n    \n    w = 2*p[-len(p)//2 : ].sum()*dx        \n    ax.plot(x, w * 0.5*np.sin(x),       c=\"darkred\",   lw=1)\n    ax.plot(x, (p - w * 0.5*np.sin(x)), c=\"darkgreen\", lw=1)\n         \n    ax.set_title(f\"angle errors: mean={np.mean(err):.3f}, median={np.median(err):.3f}, w={w:.3f}\")    \n    ax.set_xlabel(r\"$\\Delta \\Psi$ (rad)\");\n    \n#-------------------------------------------------------------------------------    \n\ndef plot_azimuth_error(err, ax,  bins = 100):    \n    \"\"\" Plot a histogram of azimuth errors and calculate the statistics and the share w of 'bad examples' \"\"\"\n    p, x, dx = plot_hist_error(err, ax, bins)    \n            \n    w  = 2*p[-len(p)//2 : ].sum()*dx                \n    v   = np.ones_like(x)/np.pi        \n    ax.plot(x, w*v,     c=\"darkred\",   lw=1)    \n    ax.plot(x, p - w*v, c=\"darkgreen\", lw=1)        \n    \n    ax.set_title(f\"azimuth errors:  mean={np.mean(err):.3f}, w={w:.3f}\")        \n    ax.set_xlabel(r\"$\\Delta \\alpha$ (rad)\");    \n\n#-------------------------------------------------------------------------------\n\ndef plot_errors(err, az_err):\n    \"\"\" Plot a histogram of angle and azimuth errors \"\"\"\n    fig, (ax1, ax2) = plt.subplots(1, 2, sharey=False, figsize=(13,4))    \n    fig.suptitle(\"Weighted Line-fit method (auxiliary==False)\")\n    plot_angle_error  (err,    ax1); ax1.set_ylabel(\"Density\"); \n    plot_azimuth_error(az_err, ax2)\n    plt.show()\n\n#-------------------------------------------------------------------------------\n\ndef angles_hist(n_pred, n_true, bins=50):\n    az_pred, ze_pred = vector2angles(n_pred)\n    az_true, ze_true = vector2angles(n_true)\n    fig, (ax1, ax2) = plt.subplots(1, 2, sharey=False, figsize=(13,4))    \n    fig.suptitle(\"Weighted Line-fit method (auxiliary==False)\")\n    ax1.hist(ze_true, bins=bins,  range=(0,np.pi), fc=\"lightgreen\", density=True, histtype='step');    \n    p,_,_=ax1.hist(ze_pred,  bins=bins, range=(0,np.pi), fc=\"lightblue\", density=True); ax1.grid(ls=':'); ax1.set_xlabel(\"zenith\")    \n    plt.legend([\"true\", \"pred\"], frameon = False, fontsize = 10)\n    \n    ax1.set_title(fr\"$\\theta$=0:  {100*p[0]/p.sum():.2f}%,   $\\theta$=pi:  {100*p[-1]/p.sum():.2f}%  (bins={bins}, $d\\theta$={np.pi/bins:.3f})\")\n    \n    ax2.hist(az_true, bins=bins, range=(0,2*np.pi), fc=\"lightgreen\", density=True, histtype='step');\n    ax2.hist(az_pred, bins=bins, range=(0,2*np.pi), fc=\"lightblue\", density=True, alpha=1); ax2.grid(ls=':'); ax2.set_xlabel(\"azimuth\")  \n    \n    plt.show()\n\n#-------------------------------------------------------------------------------\n\nplot_errors(ang_err, az_err)    \nangles_hist(-u, n)    \n","metadata":{"execution":{"iopub.status.busy":"2024-05-12T12:30:33.61234Z","iopub.execute_input":"2024-05-12T12:30:33.612939Z","iopub.status.idle":"2024-05-12T12:30:35.682148Z","shell.execute_reply.started":"2024-05-12T12:30:33.612895Z","shell.execute_reply":"2024-05-12T12:30:35.680701Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Plot some statistics","metadata":{}},{"cell_type":"code","source":"if CFG.show_stat:\n    agg_df = agg_features(pulses_df)\n    info(f\"single string events: {100*len(agg_df[agg_df.lines == 1])/len(agg_df):.3f}%\")\n    num = 11 if CFG.drop_aux else 50\n    hist(agg_df.lines, agg_df.lines[agg_df.lines <= num], pref=\"Number of strings:\", bins_sub=num)\n    hist(agg_df.doms, agg_df.doms[agg_df.doms <= 40], pref=\"Number of doms:\", bins_sub=40)\n    hist(agg_df.pulses, agg_df.pulses[agg_df.pulses <= 100], pref=\"Number of pulses: \", bins_sub=40, digits=0)\n    num = 1 if CFG.drop_aux else 3\n    hist(agg_df.t, agg_df.t[agg_df.t <= num], pref=\"Mean time of pulses: \", bins_sub=40, digits=2)\n    hist(agg_df.q, agg_df.q[agg_df.q <= 3  ], pref=\"Mean charge of pulses: \", bins_sub=40, digits=2)","metadata":{"execution":{"iopub.status.busy":"2024-05-12T12:30:35.684176Z","iopub.execute_input":"2024-05-12T12:30:35.684733Z","iopub.status.idle":"2024-05-12T12:31:51.128531Z","shell.execute_reply.started":"2024-05-12T12:30:35.68467Z","shell.execute_reply":"2024-05-12T12:31:51.126976Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"## Parameter optimization","metadata":{}},{"cell_type":"code","source":"from scipy.optimize import minimize\n\ndef fun(x):\n    pred_df = line_fit(pulses_df, lm=x[0])    \n    u = pred_df[['ux','uy','uz']].to_numpy()\n    n = true_df[['nx','ny','nz']].to_numpy()    \n    err, az_err, ze_err = angle_errors(-u,n)                     # inverse vector!\n    res = np.mean(err)\n    print(f\"err={res:.5f}  (lm={x[0]:.3f})\") # gm={x[1]:.3f}\n    return res\n\nif True:  # change this to True to run the optimization\n        x0 = [1]\n        res = minimize(fun, x0, method='Nelder-Mead', tol=1e-4)\n        print(f\"Optimized lambda value: {res.x[0]:.3f}\")\n        print(f\"Minimum mean angle error: {res.fun:.5f}\")\n","metadata":{"execution":{"iopub.status.busy":"2024-05-12T12:31:51.130787Z","iopub.execute_input":"2024-05-12T12:31:51.131996Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from sklearn.linear_model import LinearRegression\n\ntrain_data = prepare_batch(pulses, target_angles, sensors, batch_id='train')\n\n# Prepare the data\nX = train_data.drop(['target_angle'], axis=1)  # features\ny = train_data['target_angle']  # target variable\n\n# Create a linear regression model\nmodel = LinearRegression()\n\n# Train the model\nmodel.fit(X, y)\n\n# Make predictions on the test data\npred_test = model.predict(pulses)","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# from sklearn.model_selection import train_test_split\n# from sklearn.linear_model import LinearRegression\n# from sklearn.metrics import mean_squared_error\n\n# # Step 1: Prepare Data\n# # Assuming X contains the features and y contains the target variables (unit direction vectors)\n# X = agg_df[['lines', 'doms', 'core', 'pulses', 't', 'q']]  # Example features\n# y = true_df[['nx', 'ny', 'nz']]  # Example target variables\n\n# # Step 2: Split Data\n# X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)\n\n# # Step 3: Train Model\n# model = LinearRegression()\n# model.fit(X_train, y_train)\n\n# # Step 4: Predict Angles\n# y_pred = model.predict(X_test)\n\n# # Step 5: Evaluate Model\n# mse = mean_squared_error(y_test, y_pred)\n# print(\"Mean Squared Error:\", mse)\n","metadata":{"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"","metadata":{},"execution_count":null,"outputs":[]}]}