Dale Lehman writes:
Recently released preprint regarding COVID/FLU vaccines and potential risks of stroke in the aged population:
https://www.medrxiv.org/content/10.1101/2023.10.10.23296624v1.full
I haven’t read it carefully and my first impression is that it is a worthy effort and appropriately caveated. Also, the conclusions seem reasonable given the observational nature of the study and limitations of the data. My questions involve 2 things. The subgroups the examine are 65-74 years old, 75-84, and >85. I’ve seen these types of binning common in medical studies. But why? The actual ages are certainly known, so why not treat age as a continuous variable? So, of course I look for the data, thinking that I can see if treating age as continuous reveals anything of interest. Here is the link for the data:
https://www.medrxiv.org/node/714804.external-links.html
It says “All data produced in the present work are contained in the manuscript.” That seems like a worthless statement to me – it is virtually true of all papers that the data produced in the study is in the paper. The underlying data, apparently, is not addressed. No doubt it is protected for the usual privacy reasons (though I really don’t see the dangers in having this data on the anonymized 5 million+ Medicare recipients), but can’t they at least talk about the data used in the study in the data availability statement?
My reply: we discuss this general point in chapter 12 of Regression and Other Stories. Short answer is that discrete binning isn’t perfect but it’s transparent and can be better than simple linear model. One thing that people don’t always realize is that you can do binning and linear together, for example in R, y ~ z + age + age.65.74 + age75.84 + age.85.up, which can have the best of both worlds. The fitted model looks kind of goofy: it’s a step function with a slope, so has a bit of a sawtooth appearance, but gives some of the flexibility of binning along with doing something with trends within each age category.
P.S. The purpose of the above post is not to recommend fitting a model using regression with age indicators and linear age. Such a model is not even continuous; it’s a ridiculous model for any data whatsoever!
My point is that if you’re fitting a piecewise constant model (regression on age indicators), then you can add linear age to this model for free, pretty much.
It’s my impression that many users don’t realize that. Instead, they think that if they’ve “controlled for age” in their analysis, they’re done. My above post makes the general point that you can keep improving a model or an adjustment, along with the specific point that if you have age indicators you can also add linear age on top of that.
“I’ve seen these types of binning common in medical studies. But why?”
I don’t know but it’s common in lots of fields outside of medical studies, including my own field of ecology. I’ve sat in multiple grad committee meetings where a student will present continuous data, like temperature or water depth, and a faculty member will suggest that they split their data into “high” and “low” so that they can then compare things with a t-test or ANOVA. It’s maddening.
I also imagine that the study Dale brings up has a better justification than what I wrote above. Maybe they chose those brackets to match with some treatment policies? i.e., if patient is >65 do X, if >85 do Y?
Your point actually highlights one of the things that bothers me most about the treatment of age in medical studies. Many of the algorithmic guidelines have fairly poorly chosen bins. PSA is one example. The “normal” range for PSA is specified for ages 50-60, 60-70, and >70. There are huge ranges within which PSA should be expected to rise continuously. A high PSA at age 69 suddenly becomes normal when that man reaches 70. The treatment/diagnosis guidelines actually use a single PSA reading for the 60-69 age bin, despite the fact that PSA is expected to rise continuously with age. That is then compounded when studies about prostate cancer use those age bins. So, the age bins are indeed chosen based on medical algorithms, but ones that don’t make sense in this particular context.
Having reliable data on age is also a recent phenomenon. In the 1940s about 1/3 of the population was not registered at all:
https://www.usbirthcertificates.com/articles/history-birth-certificates
Before 2006 or so I don’t think there was reason to believe a 60 yr old’s reported age was much more accurate than to the decade. Binning also reduced the size of the tables getting published when that was a concern. It looks like a cultural artifact from an earlier era that requires modernization.
I believe a similar thing occurred in the COVID vaccine. Pfizer’s vaccine was approved for 12yo and older, and then quite a bit later for something like 5-12yo and then eventually 6mo to 5yr or something along those lines (I’m pulling from memory not looking it up, but I remember the 12yo cutoff quite accurately).
If you gave a 100lb 12yo boy the Pfizer vaccine it was the same dose as a 35 year old 200 lb 6ft tall man. Now, we don’t necessarily believe that mass exactly linearly relates to intensity of response, but it certainly has some relationship, and sure enough the dose for 12yo often knocked them for a loop compared to a full adult.
At the time I argued there was no reason why we should be using cutoffs, that it made much more sense to use a continuous (dimensionless) measure of dosage, and then if needed for ease of administration round off these dosages to nice round numbers and make a simple chart (perhaps 2 charts one for women and one for men)
Anyway, the problem repeated itself with the youngest age group if I remember, where they had to re-do the study after finding that the tested dose was producing too much side-effects in the 6mo age bracket or something similar.
For some reason simple cutoffs seem intuitive to non-math people. I think perhaps because you don’t have to do any mental arithmetic, it’s easy to check whether someone is more or less than 12 years old, but I still don’t find it compelling, nomogram printouts have been in use by people for decades and decades, reading values off a chart with appropriate guide-lines on it for accuracy just isn’t that hard.
I recall that we discussed this a couple years ago. I strongly agree that it made no sense to give me and my wife the same dose. Nearly every woman I know had a day of feeling bad after getting the vaccine, while nearly no man did, and I wonder how much of that is due to the fact that women were getting a larger amount per unit body size.
That said, as we discussed back then there are good reason to discretize doses. Package individual doses and there’s no need to risk contamination by drawing multiple doses from a vial, nor waste material by only providing part of what is packaged. Still, surely they could have done S-M-L-XL or similar. Make each one 1.5x the volume of the previous, so if you run out of L you can give S+M. Something like that.
Phil, the COVID vaxx was administered by drawing each dose from a vial, so they weren’t prepackaged. And that was necessary for speed of production and for storage in deep freeze (before they did experiments allowing for higher temp storage).
It’s much easier for nurses to draw the same dose every time of course… they never have to look anything up, they don’t have to ask questions that might embarrass the patient (like, how much do you weigh?), but I still argue they should have started with a continuous curve and discretized it into about 4 different doses for males and 4 different doses for females.
For adults, they’d just take the “adult dose” for females or for males, they don’t need to ask embarassing questions, it’s probably good enough. The adult female dose would be slightly lower than the adult male dose.
It’s for children between say 6 months and 16 years where they should have had 3 different doses
From a practical perspective I’d argue to do the calculation by body measurements, and then interpolate those onto the CDC standard age charts for children: https://www.cdc.gov/growthcharts/clinical_charts.htm#Set1
Then you’d ask a child how old the child is, and there’d be a very simple chart, one for males one for females, with age groups and dose in 3 different boxes or whatever. You’d ask, look up the dose, draw it, and inject.
Do the calculations based on dimensionless ratios of body measurements, but translate that into a simple to use practical chart for the nurse.
I have had a few scenarios in my career where the mix of linear and discrete terms for age make sense. So the smooth underlying plus a few big jumps at life transitions (16 for driving, 18 for going to college, 65 for retiring). I can sometimes use splines to get what I want, but it is not crazy to have the discrete jumps make sense. (I prefer splines + dummy for regression discontinuity designs.)
That makes a certain amount of sense. But I’m not seeing the sense of decades: 55-65, 65-75, 75-85, etc. What is the logic behind those bins meaning anything? Andrew’s reference to chapter 12 in Regression and Other Stories adds to my confusion. In one section, they say “We avoid discretizing continuous variables…” while the next section says “In some cases, however, it is convenient to discretize a continuous variable if a simple parametric relation does not seem appropriate.” The example is then given of using 4 indicator age ranges to capture different generational patterns. I actually find the idea of “generational patterns” somewhat bothersome, as I think the generations are artificial constructs. It seems to me that using age as a continuous variable should reveal whether there are such generational patterns. Or, if you use both, as Andrew suggests above, then you can see if generations make any sense. But to replace age as a continuous variable with a discretized version strikes me a the second step after determining that the binned pattern makes sense. And I’d want to see the first step that justifies the binning to begin with.
Oh to be clear Dale agree the age bins in the paper are silly. I think very likely to at least be linear trends (and maybe some variance in the lower end, people who select into Medicare right at 65 have different baseline health than those who wait a few years, but that in my experience washes out around early 70’s where I have looked). Figure 3 in the paper is a bit noisy, the subgroup N’s only tend to be around 500 people in those age bins.
It is a common thing in medical papers that I know Frank Harrell harps on, so I am guessing the authors do it just because it is convenient. (I mean someone may have said I need to bin for confidentiality, but I don’t see mention of that in the paper.) I think sometimes people do things not because they make sense statistically, but because they only know how to present results a certain way in a table (or here in Figure 3 with the effects per bin). I would rather they just pick one model and do various model based summaries, instead of the progressive descriptive breakdowns. (I am wondering here the days after 1-21 vs 22-42 days, I think that could maybe be a bigger binning deal than age.)
Just giving an example in general I think Andrew’s point of linear + dummies makes sense. Somewhat off-topic, I not un-commonly hear people say in regression contexts the right hand side should not be collinear, so many implicitly would think you can’t do what Andrew suggests here. Which I think is a muddled confusion between causal inference (and maybe some inverting matrix computation error stuff). But people were taught a list of “assumptions” and one point and just parrot that list afterwards.
Dale,
I wouldn’t want to argue that binning is likely to be optimal for a final model, so I agree with the general thrust of your comment. But binning can be very useful in exploratory modeling and might be all you need (although not optimal) for an entire modeling exercise.
Suppose you’re working with cancer death rate, and you start with a linear model. The probability that a given person dies of cancer in a given time interval is highly age-dependent, so you aren’t going to want to just put “age” as a linear predictor and leave it at that. You can put in “age” and “age^2”, but then you’re assuming a quadratic relationship, which might be too strong (or, perhaps, too weak). One of the simplest and most interpretable options is to create some age bins. You could just have an indicator variable for each bin, or you could fit a piecewise-continuous risk function such that the risk increases linearly within a bin, but the slope in each bin is different.
There are other things you can do, that would avoid binning — for instance, you could make up a multi-parameter risk function that increases monotonically with age, and fit the parameters — but it might not be easy to find a function that is sufficiently general that it can really fit the data, without being too sensitive to outliers or other issues. A linear model with bins is very easy to specify, calculate, and interpret, so it can be a great way to start even when the data are continuous. (This isn’t just speculation on my part, it has been my experience in a few cases, such as looking at the relationship between worker productivity and indoor carbon dioxide concentration, the latter having a direct effect on worker performance but also serving as a proxy for ventilation rate per worker).
My opinion is it’s better to know something about how to build continuous functions from some basis set. My own preference is either Chebyshev polynomials, Fourier series, Radial Basis Functions, or Splines.
(in some sense these are all closely related, Chebyshev polynomials are a transformation of a Fourier series, RBFs are a generalization of Splines, and all of them are a generalization of the notion of a linear combination of basis functions).
Often data is binned for confidentiality. For example from the US WONDER deaths data they will only give you data on number of deaths in each category when the number is greater than 10, so if you want info by age, you bin by decades and then get counts, but only if the count is at least 10.
Even with binned data, I’d rather see it represented as a function of age with 3-4 coefficients in a Chebyshev series than with 10 coefficients one for each decade, with no connection between decades at all. If death rate is high for people in their 60’s we expect it to be a nearby value for people in 70’s not plummeting to almost 0. But a coefficient per bin means there is no required relationship between the values in each bin.
Fundamentally I think statistics textbooks teach the “dummy for each category” methodology because either the authors don’t know any better or they believe the audience doesn’t understand the concept of infinite dimensional vector spaces and series representations of functions, and they don’t know how to teach it to the audience in a compact amount of time.
The other thing I should say is that something like this
y ~ a + b*x + c*x^2 + d*x^3 + e*x^4
or
y = a + b * phi1(x) + c * phi2(x) + d * phi3(x) + e * phi4(x)
where the phi1..phi4 are some chosen nonlinear functions of x,
is a *linear* model.
note that the unknowns, which are the a,b,c,d,e appear *linearly* in this model and that the equation can be written as a matrix multiplication:
y = X * (a,b,c,d,e)
where y is a column matrix of observations, X is a matrix with 1, x, x^2, x^3, x^4 as the columns, with each row being a different set of covariate measurements, and (a,b,c,d,e) is a column vector and * means matrix multiplication.
people often believe that y = a+bx is the only thing you can call a linear model because the *x* appears linearly and the graph is a line.
This is a massive misunderstanding in applied stats.
Daniel,
I would guess that most readers of this blog know that if you include age and age^2 in a linear model, it’s still a linear model. But I might be wrong!
The fact that you bring up Chebyshev, Fourier, RBF, and splines as possible families of continuous functions suggests — correctly, I think — that there is no one-size-fits-all solution, even if you think the ‘true’ function is pretty smooth. A piecewise-linear function is never right but will very often be close enough, and you can specify it in a minute. In almost all cases I would recommend it as an initial step to building a model (after all, why not?) and in many cases I’d say you can leave it there. From age 50-65, cancer risk increases at x% per year; from 65-75 it’s y%; and above 75 it’s z%. In many applications that is all you need.
I’m not saying it is _always_ all you need.
Phil,
Your piecewise linear is of course a kind of spline. So you’re saying splines, but of first order. The big problem with linear order splines is that they’re continuous functions but their derivatives aren’t continuous, and that’s kinda annoying. Like, if you want to talk about the behavior of the function in a certain range, if that range includes a knot point then the function behaves differently on one side of that knot than the other.
The usual cubical splines have continuous derivatives. Cubical splines have 4 coefficients for each section (0,1,2,3 order) which allows them to match both a given value and derivative at each end.
Chebyshev and Fourier and some RBFs are examples of *global* basis functions, and typically are infinitely smooth everywhere. That’s often a desirable property because you can then put priors on the average first, second, third derivatives etc to force them to have more or less “wiggly” behavior and soforth.
Chebyshev is for finite domains (mapped into [-1,1]), Fourier is for periodic functions on a finite domain (mapped into [0, 2pi] or [-pi,pi]) examples might include seasonality in timeseries, RBFs is for stuff on irregular domains in multiple dimensions, like a function defined on a banana shape or whatever.
So what determines which you use is more about the domain of the function and potentially how easy it is to implement in software.
This kind of topic should be about chapter 2 or 3 of “applied statistics” textbooks. It should come BEFORE probability and normal distributions etc.
Daniel,
Yes, piecewise-linear is a special case of splines. I think it’s a great default: extremely easy to implement, robust, easy to plot, and easy to interpret the coefficients. If you prefer cubic splines or something, that’s fine with me, no skin off my back! I’ve used splines many times over the years and I’m sure I will in the future too; I also often use loess and sometimes use supersmoother.
But I will continue to use piecewise-linear when I think it’s useful, which is often. Not only that, I think a lot of other people would be well served by using it too. You can tell me you don’t find it useful, and I believe you, but don’t try to tell me _I_ don’t find it useful!
Just to follow up on Phil’s comment: the purpose of my post above was not to recommend fitting a model using regression with age indicators and linear age. Such a model is not even continuous; it’s a ridiculous model for any data whatsoever!
My point is that if you’re fitting a piecewise constant model (regression on age indicators), then you can add linear age to this model for free, pretty much.
It’s my impression that many users don’t realize that. Instead, they think that if they’ve “controlled for age” in their analysis, they’re done. My above post makes the general point that you can keep improving a model or an adjustment, along with the specific point that if you have age indicators you can also add linear age on top of that.
I added the above as a P.S. just in case any new readers encounter this post.
Phil, I would never presume to tell you what you find useful!
I also agree with you that piecewise linear is not a terrible default. It’s better than piecewise constant, or than piecewise constant plus a linear component that’s discontinuous at the boundary of each region.
Still, the thing that gets me is that the general issue of representing functions in terms of a linear basis expansion has been understood well for at least a hundred years, Fourier invented the Fourier series around 1807!
Chebyshev who described Chebyshev polynomials died in 1894. Cornelius Lanczos described the importance of Chebyshev polynomial approximation theory in 1952
RBFs were studied in 1970’s and 1980’s and were being used to do map topography in the 1990’s or whatever.
(all these dates are from wikipedia articles).
The point is, well before digital computers were even available to everyday statistical data analysts (which I’d say became broadly the case in the late 70’s or early 80’s), the theory of approximation of functions by basis expansion had been worked out. Splines, Chebyshev Polynomials, and Fourier series are between 80 and 200 years old, RBFs are at least 40 years old.
So, why aren’t all of them basically accessible with something *like* (though not exactly):
y ~ Fourier(x,6)
y ~ Chebyshev(x,6)
y ~ CubicSpline(x,10)
Why do we have people still doing and even *recommending*
y ~ as.factor(theyear) + as.factor(themonth)
So, actually in Julia you can do a lot of this stuff really easily, which is one reason why I like it a lot.
FastChebInterp.jl has a function chebregression
There are almost too-many spline approximation packages in Julia. I don’t want to recommend one, but anyone who wants to do it could ask on discourse.julialang.org what the recommended of the various packages are for their purpose
and it’s easy to use ApproxFun.jl in Turing.jl to get a Bayesian fit, since you can do Fun(coefs) to construct a Chebyshev poly from coefficients, and then evaluate it at all the various x values you need.
If your book or lecture notes recommends doing dummy variables for timeseries modeling (ie. piecewise constant) then it should be thrown away and you should be laughed at and if this is common in your field you should be ridiculed as not a serious field.
y ~ as.factor(year) + as.factor(week) or whatever is a huge huge red flag
In R all you have to do is:
y ~ poly(juliandate,10)
and you’ll get a 10th order orthogonal polynomial in time, so it’s not even acceptable as a hack because of lack of computing implementation. Nor should anyone be teaching
y ~ 1 + x + x^2
instead
y ~ poly(x,2)
One can certainly specify a model with categorical and linear terms together, but to answer the question of whether the measure of association between z and y is a function of age (as in the paper), one would also need to include product terms between z and (typically) all age-related terms.
Rather than speculation, why not also share those results? I seems trivial to do once you have the ability to exclude them.
Also their smallest intervals (Fig 1) are about +/- 10%, most are larger. All-cause mortality went up ~10% during covid and this was considered a big deal. So the study was not able to detect what is otherwise considered a large effect.
which results? I think the point they’re saying is that people may have tested themselves and know they had COVID but the information never makes it into the database these people are using to do their study.
The results for the people diagnosed with covid they excluded from the study. Presumably the relationship with stroke has something to do with immune response, which should in turn be related to both infection and vaccination.
I stopped following this but if they bothered to collect such data it should be published.
Another issue is that “healthy vaccinee bias”. Ie, very ill people are less likely to be vaccinated. In this case I’d guess it would mostly be waiting for a “flare-up” of whatever illness abates before getting vaccinated. So the pre-vaccination control period isn’t quite comparable to post-vaccination.
Anyway, for a number of reasons this study doesn’t seem capable of concluding one way or the other on whether the vaccination is meaningfully increasing rate of strokes.
But that is not exceptional, it is standard for these observational studies.
If they want to do it right they need some quantitative model of what they expect to happen and then look for deviations from those predictions. That is how they are used in physics.
1. Generally, if a screened patient is excluded because they do not meet eligibility criteria, no data is collected, so there would be nothing to share.
2. Just because the confidence interval is +/- 10%, it does not follow that the study failed to detect a true effect size of 10%.
It’s not just age. I do research in a medical field, and “mild – moderate – severe mania” is omnipresent.
We do often treat age as a categorical variable and this has nothing to do with whether “the actual ages are certainly known.” Rather, assuming that age has a linear effect on the outcome of interest is often clinically unjustifiable, especially when the model has a non-identity link function. In this example, the pre-print authors use Poisson regression, so to include age as a continuous variable is to impose the unjustifiable assumption that increasing age increases the number of strokes exponentially.
That’s simply not true. You are assuming that age as a continuous variable enters entirely through a linear function of age.
q = a + b * age
for example. But there is absolutely no requirement that continuous variables be used only in linear functions.
In R something like:
glm( n ~ poly(age/100,5), family = “poisson”,data=ourdata)
now age is a number approximately between 0 and 1 and n is predicted from a 5th order numerically stable orthogonal polynomial expansion.
What I said is true for age as a “linear effect.” The reason why I assumed a linear function of age is because the title of the post is “including a regression predictor linearly.” I never made any claims about other continuous transformations of age.
Orthogonal polynomials are an interesting idea and perhaps we should consider using this more often. However, I conjecture that a major reason for why we aren’t using this more often is because the effect become less interpretable for our target audience.
Age as a factor is just a terrible way to get nonlinearity. If you have people say age 60-95 now you’ve got 35 coefficients to estimate compared to say 5-6 orthogonal polynomial coefficients. And there’s zero regulation in the factor version, your model can say that the risk for a 65 year old is 30% and for a 66 year old is .002% and for a 67 year old 1.4% etc
Interpretability is not particularly good for a factor model if you want to ask questions like how does risk increase between age 60-70 unless the model happens to fit a quite smooth function. If you have small numbers of data points at each age that won’t be the case. Like if you have 100 patients that vary over 30 years with no more than 4-5 patients in any given year and some years with nothing but 1 or even 0 patients etc.
Orthogonal polynomial fits should be the default and only move away from that if you’ve got some strong reason.
To clarify, the pre-print is looking at strokes among people aged ≥65 years. It so happens that we do have evidence that the risk of stroke increases exponentially with age in this population, so perhaps this is one clinical scenario where assuming age has a linear effect on log(outcome) is more justifiable. My point is rather that justifiable assumptions of exponentiality are rare.