<p>Extracting underlying signals from imaging data is crucial for medical applications. Medical imaging data can be contaminated by outliers or heavy-tailed noise, and the irregular domains of such data pose additional challenges. Ordinary least squares (OLS) methods are highly sensitive to outliers or heavy-tailed noise, while existing robust estimation approaches encounter difficulties with medical imaging data defined on irregular domains. To this end, we develop a novel robust estimation method by integrating M-estimation with bivariate penalized splines over triangulations. Under mild regularity conditions, we establish the <i>L</i>₂ convergence and asymptotic normality of the proposed M-estimator. Simulation studies demonstrate that the proposed method significantly outperforms OLS when the errors do not follow a normal distribution, while maintaining comparable computational efficiency. Applications to imaging data from the Alzheimer’s Disease Neuroimaging Initiative (ADNI) validate the practical value of the proposed method.</p>