Everyone says not to fit high-order polynomials, but nobody says why when you are not extrapolating

There is a real problem inside the range, and it is not overfitting in the usual sense. It is that polynomial terms are not local.

A degree-9 term is a single global function evaluated across your whole range. That means a cluster of points at the far right of your data influences the fitted shape at the far left, because both are being explained by the same coefficient on x^9. Move one outlying point at one end and the curve wobbles somewhere else entirely.

That is almost never what you want. When we fit a curve we usually believe something like "the relationship is smooth" — a local claim. Polynomials deliver "the relationship is one specific global shape", which is a much stronger claim that nobody meant to make.

The practical consequence is that the fit is unstable in a way that cross-validation on a single split will not necessarily reveal. Refit on a bootstrap resample a few times and watch the middle of the curve move.

The other half of the standard advice, which usually goes unsaid: there is a better tool for the same job, so the trade-off is not "polynomial or nothing".

Regression splines and natural splines give you the flexibility without the globalness. They are piecewise polynomials joined smoothly at knots, so a point on the right only affects the fit near the right. Natural splines go further and constrain the fit to be linear beyond the outermost knots, which kills the edge behaviour by construction.

You get roughly the same expressive power for the same number of parameters, with local influence and far more stable coefficients. Once that alternative exists, "do not use degree 9" stops being a warning about danger and becomes ordinary advice about picking the better of two tools.

If you want to keep it simpler still, a moving-window smoother will answer most "what does the shape look like" questions without committing to a parametric form at all.

Add to that: the interior misbehaviour is real, not only the edges.

The classic demonstration is fitting a high-degree polynomial through evenly spaced points on a smooth function and watching the fit oscillate hard between the points near the ends — not outside the data, between observations you actually have. Increasing the degree makes it worse rather than better, which is the opposite of the intuition that more flexibility means a better fit.

So "I only predict inside my range" does not buy the safety it sounds like it does. Inside the range still contains gaps between your observations, and gaps are where the wobble lives. If your x values are dense and even you will get away with more, but that is a property of your sampling, not of the method.

One thing that is not a real objection, since it comes up whenever this is discussed: correlated coefficients.

Raw powers of x are massively collinear, so the individual coefficients are unstable and un-interpretable — but if you only want predictions that does not matter, and using an orthogonal polynomial basis fixes the numerical side anyway. It is a genuine problem for reading the coefficients and a non-problem for fitted values.

Worth separating from the wobble, which is a real problem for fitted values.