Significance: Quantitative measurement of blood oxygen saturation (sO$_2$)with photoacoustic (PA) imaging is one of the most sought after goals ofquantitative PA imaging research due to its wide range of biomedicalapplications. Aim: A method for accurate and applicable real-time quantification of localsO$_2$ with PA imaging. Approach: We combine multiple illumination (MI) sensing with learned spectraldecoloring (LSD); training on Monte Carlo simulations of spectrally coloredabsorbed energy spectra, in order to apply the trained models to real PAmeasurements. We validate our combined MI-LSD method on a highly reliable,reproducible and easily scalable phantom model, based on copper and nickelsulfate solutions. Results: With this sulfate model we see a consistently high estimationaccuracy using MI-LSD, with median absolute estimation errors of 2.5 to 4.5percentage points. We further find fewer outliers in MI-LSD estimates comparedto LSD. Random forest regressors outperform previously reported neural networkapproaches. Conclusions: Random forest based MI-LSD is a promising method foraccurate quantitative PA oximetry imaging.