Sample size determination

A sample of 100 papers was selected to evaluate the prevalence of reported linear regression assumptions. The sample size was calculated to estimate a proportion of 0.05 (5%) with a two-sided 95% confidence interval and a margin of error of 5% [14]. It was anticipated that approximately 50% of papers would contain accessible data, enabling reproducibility assessment. This study was designed as a pilot, given the current lack of empirical evidence on the proportion of published studies that are reproducible.

Selection of papers

The initial sample of 100 research articles, excluding editorials and other non-research content, was randomly selected from PLOS ONE publications in 2019 to provide a relevant snapshot of reporting practices. This time-frame also offered an unintended advantage: it was the last full year prior to the COVID-19 pandemic, thereby avoiding potential disruptions to research [17].

The selection used the ‘searchplos’ function from the ‘rplos’ R package [18] to identify articles that met two criteria: (1) the subject area included “health,” and (2) the Materials and Methods section reported the use of “linear regression.” Eligible papers were randomly ordered, and the first 100 meeting the inclusion criteria were selected. To focus on standard applications of linear regression in health research, articles were excluded if they used clustered or random effects, alternative modelling approaches (e.g., Bayesian or non-parametric methods), or applied linear regression only for secondary purposes (e.g., calibration or pre-processing). Upon full-text review, five papers were excluded as they did not report any linear regression results, leaving 95 papers for the assessment of reporting practices.

Research questions

• What proportion of analyses and papers violated statistical assumptions (independence, linearity, homoscedasticity, and normality of residuals), including issues with outliers or multicollinearity?
• How many papers reported that the authors had checked statistical assumptions, outliers, and multicollinearity?
• What proportion of analyses remained inferentially reproducible after accounting for unmet assumptions, outliers, and multicollinearity?

Evaluating Model Assumptions and Reproducibility

Reproducibility analyses were conducted by the primary author, an accredited senior biostatistician with over 20 years of experience in statistical modelling and analysis of healthcare data. To ensure the reliability of the findings, selected results were independently reviewed by three senior biostatisticians, two of whom are accredited (AStat). Agreement on the reproduced results and interpretation was reached through discussion and consensus. The analyses followed a consistent, transparent process, documenting all steps and decisions. We refer to each paper by a unique number, the full list of identifiers and corresponding DOIs are available in the GitHub data [19].

A maximum of three linear regression analyses per paper was selected to ensure the workload remained feasible. Because models within a paper typically share the same data and analytical framework, their reproducibility was expected to be highly correlated, limiting the additional insight gained from reproducing all models. For papers with more than three models, the final or primary models of interest were included first, with any remaining models selected at random. If no primary model was identified, all three were selected at random. Models that were partially reproducible were further assessed for inferential reproducibility, using the reproduced results and available data to check model assumptions and conduct sensitivity analyses.

The assumptions of linear regression were evaluated using descriptive statistics, diagnostic plots, and formal tests, with additional checks for multicollinearity, outliers, and influential observations. A simulated example illustrating the relationship between blood pressure and age is shown in Fig 1 to demonstrate how these diagnostic procedures were applied in practice. While full implementation details are provided below and on GitHub [19], this example highlights key components of assumption checking, including inspection of scatter plots to assess linearity and potential non-linear (e.g., quadratic) relationships, examination of residual plots for evidence of non-linearity or heteroscedasticity, identification of influential observations, and assessment of whether residuals are approximately normally distributed.

Figure 1. Simulated example showing a hypothetical relationship between blood pressure and age. The top-left scatter plot shows a clear positive relationship between age and blood pressure, with the blue line showing a linear fit and the red line a GAM fit. The residual plot (top-right) indicates a quadratic pattern that is not captured by the linear model but does not show evidence of heteroscedasticity. The Cook’s distance plot (bottom-left) shows that observations 139 and 191 exert greater influence on the model, though the overall Cook’s distances remain low, suggesting that outliers are not a major concern. The Q–Q plot (bottom-right) indicates that the residuals from the linear fit are approximately normally distributed.

Visualisations of the model were generated using the ‘visreg’ package in R [20]. For linear predictors, each plot displays the relationship between the predictor (x-axis) and the outcome (y-axis), with other covariates held constant. Results for each paper are available on GitHub [19]. Each regression analysis was assessed using the steps outlined below:

Is the linear structure of the regression model appropriately specified?

• Scatterplots of standardized residuals versus predicted values and versus each predictor were used to assess the linearity of the overall model. Quadratic fits were overlaid to detect potential non-linear patterns [21].
• Residual specification tests for predictors were conducted using quadratic terms. For the overall fitted values, Tukey’s one-degree-of-freedom test was used to assess curvature and non-additivity [21].
• Scatterplots of raw dependent and predictor variables were examined with overlaid linear fits and smoothed generalized additive model (GAM) fits to explore potential non-linear relationships [22].

Do residuals display homogeneity of variance?

• The plot of residuals vs predicted values was assessed for patterns such as funnelling.
• Plots of residuals vs individual variables were used to assess patterns.
• Heteroscedasticity was assessed using the studentized Breusch–Pagan test for non constant error variance [23].

Were there outliers and influential observations?

• Cook’s distance measured the overall change in fit if the ith observation is removed. Potential influential observations are identified by
• Cook’s , where n is the number of observations [24]. A threshold of 0.5 was used to identify influential observations.
• DFFITS were used to assess how many standard deviations the fitted value for an observation changes when that observation is removed from the model [24]. A practical cut-off of 1 was used to identify observations with substantial influence.
• DFBETAS were used to assess how many standard errors a regression coefficient changes when an observation is removed from the model [24]. A practical cut-off of 1 was used to flag observations with a meaningful impact.
• Influence plots were examined to identify observations with both high leverage and large studentized residuals (typically ±2 or ±3), as these may disproportionately affect the model’s estimates. Such points appear as large, dark blue bubbles in the plots [21].
• COVRATIO measures the overall change in the covariance structure of the estimated regression coefficients, reflecting changes in their precision when each observation is removed [25]. Values close to 1 indicate little influence on the model’s precision. Values below 1 indicate inflated variances and reduced precision (wider confidence intervals), whereas values above 1 indicate deflated variances and increased precision (narrower confidence intervals). A practical cut-off between 0.9 and 1.1 was applied, with values outside this range were considered to have a meaningful impact on precision.

Are the residuals of the regression model approximately normally distributed?

• A Q–Q plot of model residuals with confidence bands assessed how residuals compared to the normal distribution. Bands indicated expected variability, and points outside the bands suggested systematic deviations [21].
• Residuals were described using descriptive statistics, including mean, standard deviation, median, minimum, maximum, skewness, and kurtosis.
• Normality was assessed using the Shapiro–Wilk and Kolmogorov–Smirnov tests [26], alongside graphical and descriptive diagnostics. Recognising that such tests are highly sensitive in large samples and may detect minor deviations, while lacking power in small samples, statistical significance alone was not taken as evidence of a meaningful violation. While formal tests contributed to the overall assessment, they were also examined to understand how often they reject minor departures from normality.

Are the residuals of the model independent of each other?

• Each paper’s study design was reviewed to assess whether independence was theoretically plausible, with particular attention to potential clustering e.g., repeated measures within individuals or observations nested within hospitals.
• Serial correlation was assessed using the Durbin–Watson test [27].

Is there collinearity among the independent variables?

• Variance Inflation Factors (VIFs) were examined to assess multicollinearity [27], with values greater than five considered potentially problematic [28].
• Changes or inflation in standard errors across regression models were monitored as indicators of multicollinearity or model instability.

Inferential reproducibility

Inferential reproducibility was assessed only for models that demonstrated sufficient computational agreement with the original results, defined as numerical agreement within pre-specified tolerance thresholds based on reported decimal precision [16]. Models classified as “Reproduced”, “Mostly reproduced”, or “Partially reproduced” were included, whereas models classified as “Not reproduced” were excluded, as meaningful evaluation of modelling assumptions requires adequate correspondence between the reproduced and original analyses. Models categorised as Partially reproduced were retained because they generally reflected the core structure of the original model, allowing assessment of whether violations of assumptions could plausibly have influenced the study conclusions. Inferential reproducibility was then evaluated as a binary outcome (Yes/No), based on whether alternative modelling approaches or assumption checks materially altered the study conclusions.

Results were presented in structured HTML reports with four tabs: Original Results, Reproduced Results, Differences, and Sensitivity Analysis, allowing readers to navigate easily and compare findings. Original results were extracted from the text and tables of the published articles, with any unreported statistics left blank to highlight gaps. The Reproduced Results tab displayed re-estimated coefficients and model fit statistics (e.g., R²), along with a summary of model assumption checks. The Differences tab was based on agreement between original and reproduced results (e.g., regression coefficients, standard errors, confidence intervals, test statistics) which were assessed using tolerance thresholds based on reported decimal place precision, classifying values as “Reproduced”, “Incorrect rounding”, or “Not reproduced”. For computational reproducibility results, see Jones et al. [19].

The Sensitivity Analysis tab of the HTML report presents inferential reproducibility assessments. Continuous dependent and independent variables from the reproduced models were standardized by centering and scaling to enable scale-free comparison of effect sizes. Standardized coefficients represent the expected change in the outcome, expressed in standard deviation units, for a one standard deviation increase in the predictor. Categorical variables were not standardized and were left in their original form (dummy variables). Standardization allowed coefficients to be compared on the same scale; however, comparing coefficients was not always an appropriate indicator of non-reproducibility, particularly when key model assumptions, such as linearity or distributional assumptions, were violated.

Differences in standardized regression coefficients were evaluated by comparing the reproduced model with sensitivity analyses assessing model assumptions. Three thresholds were applied: a percentage change of 10%, and absolute differences of 0.1 (moderate) and 0.2 (substantial). The 10% threshold was included as it was specified in the study protocol [29] and has been used in several previous studies [30, 31]. However, the percentage change can be misleading when regression coefficients are close to zero, as small absolute differences can appear disproportionately large. Because this was an exploratory project, we expanded on our protocol by also including absolute differences of 0.1 and 0.2 as more meaningful standardized measures. A difference of 0.1 can be considered similar in magnitude to a 10% change for many coefficients, while 0.2 aligns conceptually with Cohen’s d, where values below 0.2 are typically regarded as small [32]. Both percentage and absolute comparisons were explored to determine which thresholds more effectively captured meaningful differences and potential violations of assumptions.

There were five general steps followed to assess inferential reproducibility:

1. Independence design/cluster effects: If analyses were suspected of violating independence and an identification or clustering variable was available, the reproduced linear model (LM) and linear mixed model (LMM) were fit to compare the fixed effects [33]. The intraclass correlation coefficient (ICC) from the null LMM was reported [34]. Maximum-likelihood AIC values were then compared; when the LMM had a lower or similar AIC, it was retained as the primary model to preserve the data’s correlation structure for coefficient comparisons.

2. Non-linearity: Linearity was assessed by adding quadratic terms (or other non-linear terms as required). Non-linear terms were retained only when they substantially improved model fit, as indicated by reductions in AIC and BIC values and by visual inspection that showed a better representation of the observed relationship.

3. Bootstrap comparisons: Confidence intervals (95% CIs) for linear models were estimated using the bias-corrected and accelerated (BCa) bootstrap method [35]. Where BCa intervals failed to converge, percentile bootstrap intervals were used. For outcomes with non-independent observations, linear mixed models (LMMs) were fitted [33], and CIs were obtained using a parametric bootstrap based on model simulations (10,000 resamples). When heteroscedasticity was suspected, a wild bootstrap was applied in both LM and LMM settings. Bootstrap distributions were inspected for bimodality and heavy tails; where results appeared sensitive to outliers, targeted robust or sensitivity analyses were conducted.

4. Comparison of results across models: Percentage and absolute changes in estimates and confidence-interval bounds relative to the reproduced linear model were summarised using thresholds of 10% change and standardized coefficient differences of <0.10 and <0.20. Coefficient direction and statistical significance were assessed for consistency.

5. Distribution check: If inspection of model residuals or mean–variance relationships suggested that the assumed error distribution was inappropriate, sensitivity analyses were conducted using alternative distributions (e.g., Gamma for positively skewed continuous outcomes, Poisson or Negative Binomial for count data), as appropriate. Model assumptions, including the influence of outliers, were reassessed following refitting [36].

To address assumption violations, we prioritised simple, consistent approaches that enabled meaningful comparisons across models. Bootstrapping was the primary method used to account for violations of the normality assumption, with wild bootstrap [37] applied when heteroscedasticity was suspected. When linearity was questionable, models were refitted with polynomial terms, and model fit statistics such as R² and AIC/BIC, were compared rather than coefficient estimates. In cases where residuals were non-normal or showed a clear mean–variance relationship, models based on alternative theoretical distributions (e.g., Poisson, Gamma) were fitted [22]. Model fit was compared using the AIC and BIC. When multicollinearity was detected, regularisation methods such as ridge regression were applied to obtain more stable and reliable coefficient estimates [38].

When studies reported removing outliers, we did not re-include them. This decision was based on two considerations: documentation surrounding outlier removal was often poor, and data for the outlier may not have been included in the dataset. To ensure valid comparisons of model fit, models were reproduced using the same dataset used for the original analysis. Recommendations regarding outlier handling were noted when relevant. Where necessary, alternative models were compared using information criteria, interpreted according to standard guidelines. For AIC, a reduction of less than 2 indicates support of a similar model, a reduction of 4–7 indicates improved fit, and a reduction greater than 10 indicates substantially better fit [39]. For BIC, Δ BIC values of 0–2 indicate weak evidence, 2–6 positive evidence, 6–10 strong evidence, and values greater than 10 very strong evidence in favour of the model with the lower BIC [40]. These thresholds are intended as general guidelines and were interpreted in the context of model purpose, plausibility, and diagnostic assessment.

Due to our limited contextual knowledge of the original studies and our decision not to contact the authors, we retained the variables reported in each model, even when the models appeared to be overfit. In such cases, a recommendation to reduce the number of variables was noted. Models were classified as not inferentially reproducible when insufficient detail or key variables related to the study design were omitted, such as variables involving clustering, randomisation, or interactions. While efforts were made to standardize the overall approach, statistical analysis is inherently context-dependent. No single method is appropriate in all situations, and care was taken to select the most suitable approach for each case.

Statistical methods

The purpose of this study was descriptive rather than hypothesis test-driven. Descriptive statistics, including frequencies and percentages, were reported. Model-level classification of inferential reproducibility (Yes/No) was based on defined decision criteria. A standardized difference of 0.1 was used as the primary quantitative threshold, while a 10% change in key estimates and a standardized difference of 0.2 were reported descriptively for comparison. These quantitative criteria were considered alongside assessment of model assumptions (e.g., linearity and independence); the full decision framework is described above in the Inferential reproducibility section.

At the model level, up to three models were assessed for inferential reproducibility per paper, resulting in clustering of models within papers so independence could not be assumed. Accordingly, a Bayesian mixed-effects logistic regression model for the binary outcome of inferential reproducibility was fitted, with a random intercept for paper, to estimate predicted probabilities. Weakly informative priors were used to stabilise estimation given the small number of papers [41]. The intraclass correlation coefficient (ICC) and its 95% credible interval were derived from the posterior distribution of the paper-level random effect using the standard latent-scale formulation for logistic mixed models [42]. Model adequacy and convergence were assessed using standard Bayesian diagnostics, including rank-normalised convergence statistics and trace plots [43]. At the paper level, inferential reproducibility was summarised as the proportion of papers in which all assessed models met the reproducibility criterion, with 95% Wilson confidence intervals used to quantify uncertainty around this estimate. R version 4.4.2 was used for all analyses [26].

Comparison to other studies

A well-known case highlighting the consequences of inferential failure comes from economics. Reinhart and Rogoff [48] reported that countries with high debt-to-GDP ratios experienced lower economic growth, a finding that influenced global austerity policies. A spreadsheet error later uncovered by a graduate student revealed major flaws in the analysis [49]. These austerity measures led to reduced public spending in several countries, including cuts to health and social care budgets. In the UK, for example, such cuts were linked to excess deaths and widening health inequalities [50]. This case illustrates how analytical errors can shape policy with real-world health consequences, reinforcing the importance of reproducibility in evidence-based decision-making.

It is important to note that our study focused on statistical assumptions and did not consider variable selection, which can also affect reproducibility. This is clearly demonstrated by Silberzahn et al. [51], who asked 29 research teams to analyse the same dataset and research question, these types of studies are known as “Many Analyst” projects. Despite using the same data, teams applied 21 different covariate sets and a range of analytic methods, including linear regression, multilevel models, and Bayesian approaches, producing a large range in effect sizes. While 69% of teams reported statistically significant results, their conclusions varied, and choices regarding model type and handling of non-independence (e.g., variance components, clustered standard errors, or fixed effects) influenced the outcomes. The authors found that statistical expertise did not fully account for this variability, concluding that even well-intentioned analysts may reach different conclusions due to methodological flexibility in complex datasets.

Another Many Analyst project conducted by Gould et al. [52] showed how analytical decisions can substantially influence research findings, even when using the same data.

In their study, 174 analyst teams analysed two ecological datasets on blue tit nestling growth and Eucalyptus seedling recruitment. The resulting effect sizes varied widely, with some teams finding strong effects, others near zero, and some even in the opposite direction. These findings further highlight how different analytical choices can lead to markedly different conclusions.

These many-analyst studies, like our findings, highlight that reproducibility is not just about data access but also how analytical decisions, such as model structure and handling of assumptions, can substantially influence conclusions [53]. While our study focused on inferential reproducibility, others have argued that reproducibility alone is insufficient. Munafo and Smith [54] advocate for triangulation using multiple, independent approaches to answer the same question, since replication can still reinforce flawed assumptions [55]. Although debate continues about the scale of the reproducibility “crisis,” there is broad agreement that improving methodological rigour and statistical reasoning is essential [56]. Fanelli [57] proposes reframing the issue as an opportunity for researchers to take greater responsibility for research quality.

Interpreting Inferential Reproducibility: Beyond Statistical Significance

Our study contributes to the limited literature on inferential reproducibility by examining whether violations of linear regression assumptions influence study conclusions. We identified a number of papers in which linear regression was applied inappropriately, resulting in confidence intervals that were frequently too narrow. This underestimation of uncertainty suggests reduced precision and overstated confidence in the reported effects. In several cases, appropriately addressing these violations required fitting alternative models, such as replacing a linear model with a Tweedie or Inverse Gaussian model. A likely contributor to poor inferential reproducibility is limited understanding of model assumptions. Many researchers either do not formally assess assumptions or rely on statistical tests without recognising their limitations [9]. In our previous work, drawing on the larger random sample of 95 papers, we found that the most common misconception was assessing normality of the outcome rather than the residuals [14]. Of the 28 author teams who reported checking normality, only five evaluated residuals, highlighting a substantial gap between theoretical knowledge and applied practice and reflecting persistent misunderstandings about the role and interpretation of regression assumptions.

One example from the current study where a violation of normality should have prompted reconsideration of the outcome distribution is found in paper 22. While two of the three models showed only mild deviations from normality of residuals, one showed a substantial deviation. Inspection of the outcome data revealed that nearly a quarter of the values were zero. Several modelling approaches are available for handling zero-inflated data. Two types of models address structural zeros, while a third treats zeros as part of a continuous mixture [22]. For instance, imagine a study measuring post-treatment quality-of-life scores on a 0–100 scale, where zero represents no impairment. Thirty percent of patients report a score of zero, while the remaining 70% have scores between 1 and 100. This pattern fits a Tweedie distribution, a flexible family that includes the Poisson and Gamma distributions as special cases and accommodates a mean–variance relationship in which the variance increases with the mean. Alternatively, a hurdle model is appropriate when zeros arise solely from a structural mechanism, for example, when a patient did not experience the condition being assessed and therefore could not have an impairment. A zero-inflated scenario occurs when both structural zeros (patients who could not have impairment) and sampling zeros (patients who could have impairment but reported none) are present in the same dataset. The complexity of model choice and interpretation highlights the importance of consulting an experienced statistician, as papers involving statisticians tend to have fewer errors, including inappropriate test selection [58].

Assessing whether model assumptions are valid is often more similar to a mixed-methods approach than a purely quantitative process. As in qualitative research, where triangulation of multiple data sources can assess whether they tell a consistent story, assumption checking should integrate both quantitative indicators and contextual understanding. For example, knowledge of the study design is essential when assessing independence. Quantitatively, the model can be examined for the level of correlation between measurements using the intra-class correlation coefficient (ICC); the AIC can inform model complexity; and effect sizes and confidence intervals can be compared across models to assess changes. While one of these may be the dominant consideration in some cases, a full evaluation typically requires combining all sources of information to make an informed decision. For instance, when the ICC is relatively low, researchers might still choose to use a linear mixed model rather than a simple linear regression to remain consistent with the study design. However, a high ICC indicates that simple linear regression is inappropriate. We observed this in paper 16, where both models assessed had high ICCs (0.26, 0.35), indicating substantial within-cluster correlation.

In our study, p-value classifications were relatively stable, with few changes from statistically significant to non-significant and vice versa. However, this finding should be interpreted in context and cannot be directly compared to studies reporting inconsistencies in p-values [59], as our analyses were restricted to models that had already been found to be computationally reproducible. Some may argue that if statistical significance remains unchanged, there is little cause for concern. Yet statistical significance is a poor arbiter of importance. Small studies often remain non-significant despite meaningful shifts in effect magnitude or precision, whereas large studies can retain statistical significance even when changes are trivial or clinically irrelevant. Our assessment showed that many authors framed their conclusions primarily in terms of statistical significance rather than effect magnitude or precision [15]. Although significance categories were relatively stable, confidence interval widths were frequently not reproducible. This indicates that unchanged statistical significance does not imply inferential stability. Precision is particularly important when estimating clinically relevant change, as the width of the confidence interval helps define the range of plausible effects and therefore the degree of certainty when weighing potential benefits against risks for new treatments.

To evaluate changes in effect size, we examined three thresholds. A threshold of 10% change in estimates flagged many inconsequential differences, and we do not recommend this approach for assessing inferential reproducibility. In contrast, absolute thresholds of 0.1 and 0.2 for changes in standardized regression coefficients performed more appropriately. The lower threshold (0.1) was more consistent with cases in which model violations were present. For example, in one instance with an ICC of 0.26, the change in the range of the standardized confidence interval was approximately 0.15. A simulation study to evaluate the utility of these thresholds under controlled conditions would be a valuable next step.

Challenges in Assumption Checking

An important consideration in evaluating model assumptions is that the process requires judgment and visualisation, rather than relying solely on strict, hypothesis-testing-based rules, which can sometimes be misleading [60]. Even statisticians may reasonably disagree on the best approach. Multi-analyst studies have demonstrated how different analytic decisions can yield varying results [52, 51]. It is also important to note that residuals are estimates of the true errors [61]. Overly flexible approaches, such as fitting Generalized Additive Models (GAMs) to residuals can overfit noise, mistaking random variation for structure [62, 63]. While minor violations of assumptions may not always justify alternative models, we observed substantial violations in several papers in which linear regression was inappropriate yet used without question.

These challenges are compounded by the way statistics is commonly taught. Introductory courses often present methods such as t-tests, ANOVA, and linear regression in a simplified, rule-based manner, typically applied to clean, “well-behaved” data [64]. These models are rarely introduced as part of a unified framework, namely the General Linear Model (GLM), leading students to view them as unrelated procedures rather than variations of the same underlying structure [14]. Measures such as R² are frequently taught as key indicators of model quality, despite offering limited insight into model appropriateness or reliability. Broader tools for model assessment, such as the Akaike Information Criterion (AIC), may be less often introduced, leaving researchers with a limited understanding of how to evaluate fit or compare models.

Moreover, violations of statistical assumptions are often treated as technical checklist items rather than potential indicators of model misspecification. For example, apparent non-normality of residuals may reflect omitted variables or an incorrectly specified relationship between variables rather than a problem with the outcome distribution.

Similarly, outlier handling in some textbooks and basic statistics courses is often taught using simple cut-offs (e.g., three standard deviations from the mean), rather than through influence diagnostics such as Cook’s distance, DFFITS, or DFBETAS. While these measures focus on changes in regression coefficients, diagnostics that assess the precision of estimates, such as COVRATIO, which evaluates the inflation of standard errors and the potential widening of confidence intervals, are often neglected, leaving applied researchers without a complete picture of how outliers affect inference [25]. Many textbooks also provide conflicting advice on whether to exclude outliers, often without evaluating the appropriateness of removal based on influence or its impact on model validity [65].

Although excluding outliers without strong justification is generally inappropriate, influential observations should be explicitly examined in sensitivity analyses to assess their impact on estimates and inference. Approaches such as robust regression and bootstrapping are commonly used in this context. Robust regression downweights influential observations during model fitting, thereby reducing their impact on coefficient estimates and associated inference. Bootstrapping can reduce reliance on strict distributional assumptions by resampling the observed data; however, it assumes that the fitted model is correctly specified and that influential observations are representative of the underlying data-generating process. Neither approach addresses fundamental model misspecification, such as incorrectly modelling the relationship between variables.

When such influential observations were identified, their impact was evaluated to determine whether they distorted parameter estimates or reflected broader model assumption violations. In paper 46, we resolved linearity and outlier issues by log-transforming both the predictor and the outcome, but high-leverage observations continued to have a disproportionate influence on the fitted slope. Robust regression was then applied to the log–log model to down-weight these observations during estimation, producing coefficients that more closely reflected the central trend of the data. A similar issue was observed in paper 32, log-transformation similarly was used to resolved departures from linearity and substantially reduced the influence of extreme observations. However, residual heteroscedasticity persisted, indicating a violation of the constant-variance assumption. Robust standard errors were therefore applied to the log–log model to obtain valid inference without altering the point estimates. These examples highlight the distinct roles of robust-based methods: robust regression is appropriate when influential observations distort parameter estimates, whereas robust standard errors are appropriate when variance assumptions are violated but the model is otherwise adequate.

Some argue that many authors do check model assumptions but, due to space constraints in journals, do not report the details or even mention that checks were conducted [66]. To our knowledge, no empirical research has formally tested this claim. Our findings suggest that the proportion of papers clearly addressing assumptions is low. While some evidence suggests that a few authors assessed assumption violations, these violations were rarely reported in detail. For example, one paper stated that assumptions were checked and found to have only mild violations that did not substantially affect results. However, many of the papers we assessed for inferential reproducibility reported only limited, and often incorrect, assumption checks, despite major violations that substantially altered the precision of the reported effects.

Issues of assumption checking also intersect with broader concerns about modelling choices and reproducibility, where software defaults may shape analytic decisions. One example we observed is from paper 22, in which ethnicity was incorrectly modelled as a continuous variable. In SPSS, unless categorical variables are explicitly dummy coded, they are treated as continuous in linear regression models. Similar risks arise when stepwise regression is applied to dummy-coded categorical variables and interaction terms, where software may treat individual terms as independent. A possible example is from paper 38, where substantial effort was devoted to testing assumptions of their interaction model; yet backward stepwise regression was subsequently used, resulting in the elimination of main effects involved in the interaction and components of categorical variables. SPSS also includes a function called Univariate General Linear Model, which performs linear regression within the broader General Linear Model framework. This function automatically dummy-codes categorical variables, which researchers may not realise given its name and setup, although it does not perform variable selection. These examples highlight a key point: understanding your software’s default settings is essential. Seemingly minor choices, such as drop-down selections or default parameterisations, can substantially change the interpretation of results. Without adequate statistical training, researchers risk producing results that are computationally reproducible but inferentially misleading.

Our findings highlight the need for better statistical education that reflects the complexity and judgment involved in modelling. Qualitative evidence suggests that researchers often interpret regression coefficients through context dependent reasoning rather than as purely mechanical statistical outputs [67], reinforcing the importance of conceptual understanding over rule based application. Some authors in our sample made genuine efforts to satisfy regression assumptions, but this sometimes led to models that were difficult to interpret and overemphasised statistical significance. For example, transformations were sometimes applied and not back-transformed, limiting interpretability on the original scale. In other cases, transformations appeared to be motivated primarily by a desire to “pass” normality tests, even when deviations were unlikely to be practically meaningful. Teaching statistics within a broader modelling framework, one that emphasises critical thinking over rigid rules, would better prepare researchers to apply methods appropriately, consider alternatives such as bootstrapping or other distributions, and prioritise interpretability and practical relevance over arbitrary thresholds.

General recommendations for improving inferential reproducibility

• Consider the underlying data-generating mechanism (e.g. continuous, count, or proportions) when selecting analytical approaches. Provide exploratory visualisations in the Supplementary Information, such as scatterplots, boxplots, or other appropriate plots, to illustrate key relationships between variables included in the analyses.
• Clearly describe the study design, explicitly stating whether the study is observational or experimental, whether randomisation was used, and whether there is temporal ordering, such as in a baseline and follow-up design.
• Provide a unique observation identifier (ID) in the dataset to ensure each individual is uniquely and consistently identified.
• Consider whether the study design includes clustering (e.g. by site, school, hospital, or family). If clustering is present, include the relevant clustering variables in the dataset.
• Clearly describe and justify the modelling strategy, including the rationale for inclusion of confounders, mediators, and interaction terms, and whether a formal model-building framework was used.
• If a formal modelling process was used, the variable selection and/or removal procedures should be described in detail, including all decision thresholds and criteria. Authors should also report the number and type of candidate variables initially considered for inclusion, rather than only listing those retained in the final model.
• Evaluate the assumptions of linear regression models by examining residual diagnostics, identifying influential observations, and assessing multicollinearity among predictors. Where assumptions are violated, apply and report appropriate remedial methods.
• If model assumptions are violated, consider fitting alternative models that better reflect the underlying data structure (e.g. non-linear specifications or models with different distributional assumptions). Use appropriate model-fit statistics, such as R2, Akaike Information Criterion (AIC), or Bayesian Information Criterion (BIC), to compare competing models and justify the final model choice.
• Assess heteroscedasticity explicitly and, where present, consider using robust standard errors or alternative modelling strategies to ensure valid statistical inference. When systematic variance differences are evident (e.g., differing variances across groups or increasing variability with a predictor), weighted least squares may be considered with appropriate justification.
• Identify influential observations using appropriate influence diagnostics (e.g. Cook’s distance, DFBETAs, DFFITS, and covariance ratios), and evaluate their impact on coefficient estimates and statistical inference. Where influential observations materially affect results, conduct sensitivity analyses, such as applying robust regression methods to reduce their influence.
• Transformations may be used to address violations of linearity or homoscedasticity, but should be applied only where they meaningfully improve model adequacy while retaining interpretability. For log-transformed models, back-transformation is often required to present results on the original scale, such as geometric means or ratios (interpreted as percentage or fold changes). Researchers should consider whether a transformation is necessary and whether the resulting parameter estimates correspond to the effect of interest, recognising that transformations (e.g. log scale) change interpretation from additive differences to multiplicative or proportional effects.
• Bootstrapping can provide standard errors and confidence intervals without relying on the assumption that residuals are normally distributed. However, it does not address model misspecification, such as non-linearity. When bootstrapping is used, the resampling scheme (e.g., residual, case, wild, or cluster bootstrap), the number of resamples, and the model coefficients and confidence intervals should be clearly reported.
• Report the details of assumption checking in the Supplementary Information, including diagnostic plots and any sensitivity analyses undertaken.
• Provide a clear description of the extent and patterns of missing data, including the assumed missing data mechanism (e.g., missing completely at random, missing at random, or missing not at random). Authors should detail the analytical approach used to address missingness (e.g., complete case analysis, single or multiple imputation, inverse probability weighting), along with any sensitivity analyses conducted to assess the robustness of results.
• Share code, underlying data, and a comprehensive data dictionary with a well-documented workflow to clarify modelling decisions and support inferential reproducibility.

Limitations

While every effort was made to faithfully reproduce published analyses using the methods described in the original articles, evaluating statistical assumptions is inherently subjective. Contextual knowledge, often unavailable to external analysts, can influence judgments about the severity and relevance of assumption violations.

Consequently, reasonable differences in interpretation and analytic decisions may arise among statisticians and researchers. In addition, implicit methodological details, such as data-handling procedures or variable coding, were often not fully reported, which may affect reproducibility. Nonetheless, we maintain that research articles should provide sufficient detail for any qualified statistician to independently reproduce the analyses.

Our results should be interpreted with caution, as we only assessed inferential reproducibility in papers that were computationally reproducible, a subset that may disproportionately represent methodologically skilled authors. The primary author evaluated inferential reproducibility, with input from co-authors as required. Although our protocol initially specified a 10% change in coefficients as a threshold, this was often not meaningful when coefficients were close to zero. To address this, additional comparisons were made using cut-offs of 0.1 and 0.2 for standardized regression coefficients.

As this was an exploratory study, the methodology was refined as challenges were identified; the original study protocol is available on GitHub [29]. Overall, our findings align with previous research documenting concerns about inferential reproducibility in scientific research [3] and with limited or inconsistent reporting of model assumptions [68].