Definition of Zernike Polynomials

Zernike polynomials are a set of orthogonal functions defined over the unit disk, typically expressed in polar coordinates (ρ, θ). They form a complete basis for representing smooth, continuous functions on a circular aperture, making them ideal for modeling optical wavefronts. Each polynomial is indexed by radial degree n and azimuthal frequency m, with forms like piston (n=0, m=0), tilt, defocus, astigmatism, coma, and higher modes. Their orthogonality ensures efficient decomposition without redundancy, which is why they’re a staple in wavefront analysis.

They are useful for modeling complex aberration in the light path in the eye. Let’s summarize them by distinguishing 2 main types : Low Order Aberrations and High Order Aberrations

What Are High-Order Aberrations? Aberrations are classified into low-order (LOAs) and high-order (HOAs). LOAs include common refractive errors like myopia (nearsightedness), hyperopia (farsightedness), and astigmatism, which stem from basic shape issues in the eye and can typically be corrected with glasses, contacts, or standard surgery. In contrast, HOAs are more subtle and complex distortions, such as coma, trefoil, or spherical aberration, resulting from minute corneal irregularities, scars, or trauma. These cannot be fully addressed by traditional glasses and often lead to symptoms like glare, halos, or reduced contrast sensitivity, especially in low light.

How High-Order Aberrations Affect Light Entry in Diseases Like Keratoconus In conditions like keratoconus, the cornea thins and bulges into a cone shape, severely disrupting light entry. This causes pronounced HOAs, where light rays don’t converge properly on the retina, resulting in blurred, distorted vision. Patients experience halos around lights, glare, light sensitivity, night vision problems, and even headaches or eye strain. Daily life is profoundly impacted: driving becomes hazardous due to poor night vision, reading or using screens causes fatigue. Simple tasks like recognizing faces or navigating spaces turn challenging, affecting work, social interactions, and independence.

Transforming Raw Data from the MS39 Device into Zernike Polynomials The MS39 is a state-of-the-art corneal topographer and tomographer from CSO, combining Placido disk topography for surface mapping with spectral-domain optical coherence tomography (SD-OCT) for cross-sectional imaging. It provides detailed data on pachymetry (thickness), elevation, curvature, and dioptric power across anterior, posterior, and stromal surfaces. My pipeline processes raw CSV exports from the MS39, extracting diameter zones (2-7mm), fitting a best-fit sphere for centering, computing Zernike coefficients up to order 6, applying optical path difference corrections, and generating visualizations. This transforms noisy, polar-coordinate data into standardized, physiologically meaningful polynomials. Goal: Understanding and Predicting Complex High-Order Aberrations with Computer Vision The ultimate aim is to leverage computer vision (CV) techniques to not just describe but predict HOAs from MS39 data. By extracting Zernike coefficients as features, we can feed them into ML models to forecast disease progression, classify conditions like keratoconus, or optimize custom treatments. This pipeline serves as a feature engineering foundation, enabling CV algorithms to analyze topography maps for patterns invisible to the naked eye. Methods for Zernike Pipeline in MS39 Data Processing This section outlines the core methods of a robust pipeline for processing MS39 corneal topography data into Zernike polynomials, designed to generate ML-ready features for predicting high-order aberrations (HOAs). Each method is critical for transforming raw data into physiologically meaningful coefficients, with concise code snippets provided.

MS39 Data Processing Pipeline

1. Data Ingestion and Preprocessing

Purpose: Extracts elevation data (anterior, posterior, stromal surfaces) from MS39 CSV files across 2-7mm diameter zones, handling noise and missing values.
Why It Matters: Creates a structured dataset for ML feature extraction, demonstrating robust data handling for medical imaging. This part of the code is tailored for MS39 exports, which are polar coordinate data.

def process_single_csv_zones(p):
    r = {'filename': os.path.basename(p), 'surfaces': {}}
    for sn, sr, nr, mn, mx in segments:
        ss, m = sn.replace('elevation_', ''), lire_segment(p, sr, nr)
        if m.empty: continue
        m = m.mask((m < mn) | (m > mx)).replace(-1000, np.nan)
        r['surfaces'][ss] = {}
        for dia in [2, 3, 4, 5, 6, 7]:
            zd = extract_diameter_zones(m, dia)
            if zd is not None: r['surfaces'][ss][f'{dia}mm'] = {'data': zd}
    return r

2. Best Fit Sphere (BFS) Centering

Purpose: Fits a sphere to elevation data using least-squares optimization, centers coordinates on the sphere’s vertex, and computes residual maps. Why It Matters: This is a standard step in corneal topography analysis, ensuring that the residual maps are centered on the sphere’s vertex. Without it, the maps would only showcase standard light aberrations (like coma and spherical aberration) without accounting for the geometry of the cornea. With it, we can isolate the aberrations needed for clinical diagnosis.

def calculate_bfs(zd,rs=0.2):
    x,y,z=polar_to_cartesian_coords(zd,rs)
    if len(x)<4:return None,None
    zr=np.max(z)-np.min(z)
    er=np.clip((np.max(np.sqrt(x**2+y**2))**2)/(2*zr),6.,15.) if zr>1e-3 else 8.
    ig=[0.,0.,np.max(z)+er,er]
    res=least_squares(bfs_residuals,ig,args=(x,y,z),method='lm',max_nfev=1000)
    if not res.success or not (abs(res.x[0])<5 and abs(res.x[1])<5 and 4.<res.x[3]<20.):return None,None
    return {'center_x':res.x[0],'center_y':res.x[1],'center_z':res.x[2],'radius':res.x[3],'rms_residual':np.sqrt(np.mean(res.fun**2)) },{'n_points':len(x),'convergence':res.success}

3. Zernike Coefficient Fitting

Purpose: Fits Zernike polynomials (up to order 6) to residual maps, producing coefficients for aberrations like coma and spherical aberration. Why It Matters: Generates compact, interpretable features for ML tasks like disease classification or progression prediction.

def fit_zernike_coefficients(x,y,z,mo,pr):
    vm=np.isfinite(x)&np.isfinite(y)&np.isfinite(z)
    xv,yv,zv=x[vm],y[vm],z[vm]
    if len(xv)<10:return None,None,None
    rho,theta=np.sqrt(xv**2+yv**2)/pr,np.arctan2(yv,xv)
    i=rho<=1.001
    rho,theta,zv=rho[i],theta[i],zv[i]
    if len(rho)<10:return None,None,None
    modes=generate_zernike_modes(mo)
    A=np.vstack([zernike_polynomial(n,m,rho,theta) for j,n,m,name in modes]).T
    c,_,_,_=np.linalg.lstsq(A,zv,rcond=None)
    if not np.all(np.isfinite(c)):return None,None,None
    zf=A@c
    ve=100*(1-np.sum((zv-zf)**2)/np.sum((zv-np.mean(zv))**2)) if np.var(zv)>1e-12 else 100.
    return c,modes,{'rms_residual':np.sqrt(np.mean((zv-zf)**2)),'n_points_fitted':len(rho),'variance_explained_pct':ve}

4. Optical Path Difference (OPD) Correction

Purpose: Converts geometric coefficients to OPD units using a correction factor (-0.376) and zeros the piston term. Why It Matters: This reproduces the industry standard where the (0,0) is NOT at the center of the cornea but at the Apex (most elevated point) of the cornea

def apply_opd_and_convention_correction(coefficients_mm, modes):
    correction_factor = -0.376
    corrected_coefficients_mm = np.array(coefficients_mm) * correction_factor
    if len(corrected_coefficients_mm) > 0:
        corrected_coefficients_mm[0] = 0.0
    return corrected_coefficients_mm

5. Visualization and Output

Purpose: Saves Zernike coefficients to CSV and generates aberration maps for validation. Why It Matters: Ensures interpretability and trust in ML pipelines, crucial for healthcare applications.

def save_results_to_csv(ar, of):
    sf = {k: v for k, v in ar['files'].items() if v.get('status') == 'success'}
    if not sf: return
    all_geom_results, all_opd_results = [], []
    rm = generate_zernike_modes()
    ts = datetime.now().strftime("%Y-%m-%d_%H-%M-%S")
    for fn, fd in sf.items():
        geom_file_row = {'filename': fn}
        opd_file_row = {'filename': fn}
        for sn, sd in fd.get('zernike_results', {}).get('surfaces', {}).items():
            for zn, zd in sd.items():
                p = f"{sn}_{zn}"
                c_geom = zd.get('coefficients_geometric_mm')
                c_opd = zd.get('coefficients_opd_corrected_mm')
                for j, n, m, name in rm:
                    geom_file_row[f'{p}_zernike_z{j}_{name}_um'] = c_geom[j] * 1000 if c_geom is not None and j < len(c_geom) else np.nan
                    opd_file_row[f'{p}_zernike_z{j}_{name}_um'] = c_opd[j] * 1000 if c_opd is not None and j < len(c_opd) else np.nan
        all_geom_results.append(geom_file_row)
        all_opd_results.append(opd_file_row)
    df_geom = pd.DataFrame(all_geom_results)
    df_opd = pd.DataFrame(all_opd_results)
    df_geom.to_csv(os.path.join(of, f"batch_results_geom_{ts}.csv"), index=False, sep=';', float_format='%.6f')
    df_opd.to_csv(os.path.join(of, f"batch_results_opd_{ts}.csv"), index=False, sep=';', float_format='%.6f')

Relevance to ML: This pipeline transforms complex MS39 data into Zernike coefficients, ideal features for ML tasks like classifying keratoconus or predicting HOA progression. It showcases skills in data preprocessing, optimization, and feature engineering, directly applicable to Machine Learning Engineering.