Optimized MAP solution is wrong after clipping residuals

I've noticed some strange behavior when I'm fitting secondary eclipse light curves with noninformative priors. Like in the TESS case study, for each light curve I optimize for the model parameters, mask residual outliers, and rerun the optimizer, using the final result as the starting point for MCMC sampling. For some of my Kepler quarter light curves, interestingly, after I mask out residuals and rerun the optimizer, the MAP solution is completely wrong, while the MAP solution before clipping was a good fit. I checked that the outliers that were being removed were actual outliers (e.g., 5 out of the 4134 data points).

To fix this and have the second optimization work, I had to set its starting point to be the MAP solution from the first optimization. Is there a reason why this is necessary? It seems very unintuitive to me why the fit would get worse using the same starting point and initial priors after removing outliers. I only observed this happen in some of my light curves.

I've attached an example plot below. The GaussianProcess + light curve is plotted in orange/red.
Phase folded plot of first fit using noninformative priors:

Phase folded plot of second fit after clipping 5 residuals and using the same test point and noninformative priors as first:

Thanks in advance, and I'd be happy to elaborate!

Start by reproducing the secondary-eclipse fitting workflow described in the issue: optimize, mask residual outliers, and optimize again with the same starting point and priors. Compare the two MAP results on the attached-style case; done means the cause of the divergent fit is identified and the clipped-data optimization behaves as intended.