Purpose:Recent clonogenic cell survival and γH2AX studies suggest proton relative biological effectiveness (RBE) may be a non‐linear function of linear energy transfer (LET) in the distal edge of the Bragg peak and beyond. We sought to develop a multiscale model to account for non‐linear response phenomena to aid in the optimization of intensity‐modulated proton therapy.Methods:The model is based on first‐principle simulations of proton track structures, including secondary ions, and an analytical derivation of the dependence on particle LET of the linear‐quadratic (LQ) model parameters α and β. The derived formulas are an extension of the microdosimetric kinetic (MK) model that captures dissipative track structures and non‐Poissonian distribution of DNA damage at the distal edge of the Bragg peak and beyond. Monte Carlo simulations were performed to confirm the non‐linear dose‐response characteristics arising from the non‐Poisson distribution of initial DNA damage.Results:In contrast to low LET segments of the proton depth dose, from the beam entrance to the Bragg peak, strong deviations from non‐dissipative track structures and Poisson distribution in the ionization events in the Bragg peak distal edge govern the non‐linear cell response and result in the transformation α=(1+c_1 L) α_x+2(c_0 L+c_2 L^2)(1+c_1 L) β_x and β=(1+c_1 L)^2 β_x. Here L is the charged particle LET, and c_0,c_1, and c_2 are functions of microscopic parameters and can be served as fitting parameters to the cell‐survival data. In the low LET limit c_1, and c_2 are negligible hence the linear model proposed and used by Wilkins‐Oelfke for the proton treatment planning system can be retrieved. The present model fits well the recent clonogenic survival data measured recently in our group in MDACC.Conclusion:The present hybrid method provides higher accuracy in calculating the RBE‐weighted dose in the target and normal tissues.