What the Stars Taught Me About Fitting a Model

I am not an astronomer, but I find the science fascinating. One pleasure of building RainbowStats is that a public dataset can take me somewhere unfamiliar. I can draw a picture, ask a question, and discover how much I still have to learn.

This time it was stars. The dslabs stars dataset [1] pairs temperature with absolute magnitude. The relationship plainly bends. That raised a practical question: why should exploring an equation require more work than writing it down?

We built CROSS_SECTION_MODEL to make that step easier. Supply a cross-sectional dataset, write an equation using its column names, and estimate the coefficients. Then inspect the fitted curve, the errors, and the observations that resist the explanation.

First, understand the axes

The compilation contains 96 rows in our analysis, with star names, temperature in kelvin, absolute magnitude, and spectral type. Absolute magnitude standardizes brightness to a common distance; smaller numbers indicate brighter stars. Some star names repeat, so 96 rows should not be read as 96 independent stars.

Download the 96 plotted temperature–magnitude pairs or obtain the original CSV with names and classifications. The included pairs were exported from the chart output and retain its displayed precision.

A straight line, a parabola, and a logarithm

We fit three equations to the same observations. With m representing magnitude and T temperature, they are m = a + bT, m = a + bT + cT², and m = a + b log₁₀(T).

These are response-versus-temperature plots. They differ from an actual-versus-predicted diagnostic, whose diagonal marks perfect prediction. Their axis directions also differ from the conventional Hertzsprung–Russell diagram.

What improved?

RainbowStats results, November 7, 2026. All models use all 96 rows and the same response scale. RMSE is in magnitude units. Lower RMSE and AIC are preferable; AIC includes the fitted error variance in its parameter count.

The logarithmic equation does best of these three. Its leave-one-out error is about 12% lower than the straight line’s, and it uses one fewer coefficient than the quadratic. The estimated equation is:

That is an improvement in description, not an astronomical discovery. The R² of about 0.54 leaves substantial variation unexplained. And leave-one-row-out validation can be optimistic when another observation of the same star remains in the training sample. A stronger analysis would reconcile identities and hold out whole stars.

The missing information is interesting

Stellar luminosity depends on size as well as temperature. The physical relationship is L = 4πR²σT⁴: a star’s emitting area multiplies its radiative flux. Stars with similar temperatures can therefore have very different luminosities. The familiar main sequence, giants, and white dwarfs occupy different regions of the H–R diagram [2].

Our temperature-only equation cannot account for that structure. Nor does its fitted slope directly test the fourth-power law: that law concerns total, bolometric luminosity, while the dataset’s magnitude is not documented here as a bolometric measure. The logarithmic model is a useful descriptive candidate, not a complete physical model.

Inspect the logarithmic model’s residual plot

From stars to streams, crops, and economies

This is the part I find most useful: the tool is general. A cross section might contain stars, stream sites, farms, households, firms, or countries. Multiple explanatory columns can enter one equation, including powers, logarithms, and interactions.

For example, a hypothetical ecological dataset with the appropriate numeric columns could use:

EcoModel=CROSS_SECTION_MODEL(Eco,"column(biomass)=a+b*column(temp)+c*column(nutrients)+d*column(temp)*column(nutrients)")

The interaction asks whether the temperature relationship changes with nutrient levels. This is an illustrative equation, not an ecological finding. The outcome and sampling design still determine whether least squares is appropriate; a count or binary outcome may call for another model.

CROSS_SECTION_MODEL currently fits equations that are linear in their unknown coefficients. It can fit a+b*LOG10(column(temp)) or a polynomial, but it is not an unrestricted nonlinear parameter optimizer. Curvature in a predictor and nonlinearity in an unknown coefficient are different things.

My goal is to make the distance from a question to an examinable model shorter. The curve is a beginning. The observations it misses often suggest the next question.

Run the analysis

Stars=CROSS_SECTION("https://cran.r-universe.dev/dslabs/data/stars/csv")
M=CROSS_SECTION_MODEL(Stars,"column(magnitude)=a+b*LOG10(column(temp))")
SLIDESHOW(M)

Open the complete studyDownload the full annotated script

If the embedded app is unavailable or cramped, use “Open the complete study.” The charts and results above remain readable independently.

Appendix: sources and related studies

1. dslabs: Physical Properties of Stars. Dataset documentation and column definitions; provenance is a compilation of open-access VizieR references. Source CSV.
2. Todd Thompson, Ohio State: Stars & the Hertzsprung–Russell Diagram. Background on temperature, radius, luminosity, and stellar populations.
3. robustbase: starsCYG. A related teaching example containing 47 observations, documented as 43 main-sequence stars and four giants. Source: Rousseeuw and Leroy (1987), Robust Regression and Outlier Detection, p. 27, table 3, using data from C. Doom. It illustrates why observations far from a fitted band may belong to another physical population.
4. Gaia Collaboration (2018), Gaia Data Release 2: Observational Hertzsprung–Russell diagrams, Astronomy & Astrophysics 616, A10. Open-access record. Further reading on the much richer structure visible in large stellar samples; this post does not replicate that study.

Reproducibility notes

Charts and statistics were obtained from the RainbowStats script accompanying this post, using all 96 complete rows. No outliers were removed. The teaching compilation is not a random sample of stars. Curve comparison was exploratory; reported validation scores are not an independent final test after model selection. The fitted lines stay within the observed temperature range. Names and classifications were not used as predictors.