The present invention provides a system and method for generating intensity modulated proton therapy plans. The method includes receiving a representation of a subject that includes information related to target and non-target volumes of interest, and computing, using the representation and a proton beam configuration that describes a number of beamlets and their respective arrangement, a dose distribution for each one of a plurality of uncertainty scenarios. The method also includes computing a robustness index for respective locations in the target volumes of interest using the dose distributions, and optimizing an objective function to determine a set of intensity weights for each beamlet in the proton beam configuration, wherein the objective function is based at least in part on deviations of the robustness indices from a prescribed robustness index. The method further includes generating a treatment plan using the determined intensity weights.
