ABSTRACT We present a method to self-consistently propagate stellar-mass $\hbox{$\hbox{${\rm M}$}_{\star }$}=\log (\hbox{${\rm M}$}/\hbox{${\rm M}_{\odot }$})$ and star-formation-rate $\hbox{${\Psi }$}=\log (\hbox{${\psi }$}/\hbox{${\rm M}_{\odot }$}\, {\rm yr}^{-1}$) uncertainties on to intercept (α), slope (β), and intrinsic-scatter (σ) estimates for a simple model of the main sequence of star-forming galaxies, where $\hbox{${\Psi }$}= \alpha + \beta \hbox{$\hbox{${\rm M}$}_{\star }$}+ \mathcal {N}(0,\sigma)$. To test this method and compare it with other published methods, we construct mock photometric samples of galaxies at z ∼ 5 based on idealized models combined with broad- and medium-band filters at wavelengths 0.8–5 μm. Adopting simple Ψ estimates based on dust-corrected ultraviolet luminosity can underestimate σ. We find that broad-band fluxes alone cannot constrain the contribution from emission lines, implying that strong priors on the emission-line contribution are required if no medium-band constraints are available. Therefore, at high redshifts, where emission lines contribute a higher fraction of the broad-band flux, photometric fitting is sensitive to Ψ variations on short (∼10 Myr) time-scales. Priors on age imposed with a constant (or rising) star formation history (SFH) do not allow one to investigate a possible dependence of σ on $\hbox{${\rm M}$}_{\star }$ at high redshifts. Delayed exponential SFHs have less constrained priors, but do not account for Ψ variations on short time-scales, a problem if σ increases due to stochasticity of star formation. A simple SFH with current star formation decoupled from the previous history is appropriate. We show that, for simple exposure-time calculations assuming point sources, with low levels of dust, we should be able to obtain unbiased estimates of the main sequence down to $\mathrm{ log}(\hbox{${\rm M}$}/\hbox{${\rm M}_{\odot }$})\sim 8$ at z ∼ 5 with the James Webb Space Telescope while allowing for stochasticity of star formation.