Zero-inflated models for the evaluation of colorectal polyps in colon cancer screening studies—a value-based biostatistics practice

Introduction

Colorectal cancer (CRC) is the second leading cause of cancer-related mortality and the third most common cancer in the world. The estimated number of new cases of CRC is highest in the United States (Xi & Xu, 2021). CRC more often begins with a polyp, a noncancerous progression that develops in the large intestine. Some types of polyps have the potential to develop CRC if left untreated (Helsingen & Kalager, 2022). The leading precursor lesion of CRC is adenomatous polyps (Klos & Dharmarajan, 2016). Although CRC is a major health concern, declining trends have been observed in the incidence and mortality rates in the United States due to the increased use of colonoscopies for CRC screening and early treatment of polyps (Kusnik et al., 2023). Continuous advances in colonoscopy procedures are being made to improve the diagnosis, early detection, and treatment of polyps among at-risk populations (Pamudurthy, Lodhia & Konda, 2020). The performance matrix of colonoscopy studies typically includes the detection of polyps or adenoma counts (Amano et al., 2018; Delavari et al., 2015; Winawer et al., 2006). Unfortunately, these outcomes in colonoscopy studies are analyzed insufficiently and potentially inaccurately.

Since researchers are often interested in evaluating polyp detection rate (PDR) in colonoscopy studies, the number of polyps is analyzed in their dichotomized form, which may produce insufficient findings (Liu et al., 2020; Rutter et al., 2020). In statistical parlance, analyzing data deviated from the original form of data collection typically causes information loss and produces inefficient results (Altman & Royston, 2006; Dwivedi, 2022). Many research studies are available for evaluating colorectal polyps in colonoscopic studies. However, most of these studies analyzed colorectal polyps in their dichotomized form (Nazarian et al., 2021; Pan et al., 2020). A study highlighted the use of count data models for the analysis of the number of adenomatous polyps in a phase III randomized trial. This study demonstrated how critical it is to analyze the number of adenoma polyps instead of their categorized form (Xie & Aickin, 1997). A study established that the count of polyps follows a Poisson distribution (Emerson et al., 1993). Recently, a few studies applied the Poisson regression (PR) model to analyze the number of polyps in colonoscopic studies (Drew et al., 2016; Liu et al., 2020). However, a large proportion of colonoscopic data yields no polyps, leading to an overdispersion issue, i.e., observed variance is greater than the estimated variance from the model in analyzing polyp data (Davies et al., 2023; Hull et al., 2019; Singh et al., 2017). The overdispersion issue in the standard PR model can be addressed by using a robust variance estimation approach in the PR model using a quasi-Poisson (QP), robust Poisson (RP), overdispersion Poisson, negative binomial (NB), or generalized Poisson (GP) model instead of a standard Poisson model (Harris, Yang & Hardin, 2012; Payne et al., 2018). However, due to an excess amount of zeros (no polyps), these single-distribution models may not be reliable enough to produce accurate findings (Campbell, 2021). Although clinicians are often interested in PDR estimation and comparison in the primary analysis of colonoscopic studies, which cannot be directly obtained by simple count data models, these can be estimated indirectly via computing marginal effect sizes, which are typically useful in clinical trial studies.

While evaluating the count data between different groups, analysts sometimes use the logistic regression (LR) model by categorizing data into a binary form. However, this approach does not demonstrate whether groups are significantly different in detecting the extent of counts. Generally, analysts apply standard count models such as PR, RP, GP, or NB models to evaluate count differences between groups. However, this approach produces a potentially inaccurate effect size as well as inadequate conclusions in the presence of excess zeros (Fernandez & Vatcheva, 2022; Slymen et al., 2006). Therefore, most trials evaluating polyp count data using these methods might have yielded potentially inaccurate conclusions or inaccurate effect sizes, affecting scientific progress and planning of future trials. As per value-based biostatistics practice, multiple questions should be addressed simultaneously if the data structure allows (Dwivedi, 2022). Accordingly, we propose to analyze polyp data using mixed distribution models such as zero-inflated or zero hurdle count models, which would not only account for overdispersion due to excess zeros but also facilitate clinical interpretation of PDR along with the estimation and comparison of the extent of detected polyps. To the best of our knowledge, as of December 12, 2024, a quick search using the terms “zero-inflated” or “zero hurdle” on Google Scholar yielded 17,900 results, while these searches on PubMed yielded 2,565 results. Although inflated and hurdle methods are progressively used in medical research (Haslett et al., 2022), these methods have not been used for analyzing studies with polyp count data. There are multiple types of inflated or hurdle models, including zero-inflated Poisson (ZIP), zero-inflated robust Poisson (ZIRP), zero-inflated negative binomial (ZINB), zero-inflated generalized Poisson (ZIGP), zero hurdle Poisson (ZHP), and zero hurdle negative binomial model (ZHNB). The choice of hurdle and inflated models depends on the type of excess zeros, while the choice of specific distribution in the inflated or hurdle models depends on the over or under, or equidispersion outcome data and goodness of the fit of a specific model (Dwivedi et al., 2010, 2014b). After eliminating excess zeros from the response, if the variance of the count data is equal to its mean, then it is called equidispersion, otherwise overdispersion (mean < variance) or underdispersion (mean > variance). We demonstrate the significance of using inflated models in analyzing polyp data using our three published studies (Avalos et al., 2020; Elhanafi et al., 2017; Othman et al., 2017) as well as the step-by-step approach for selecting an appropriate model.

Materials and Methods

In this report, we used published studies to demonstrate the utility of inflated models and the interpretation of data in analyzing polyp response yielded from either observational or experimental studies. Two of the three example datasets in this current study were generated from a randomized clinical trial.

Statistical analysis

We followed the steps for selecting an appropriate model as specified in Fig. 1. Accordingly, we assessed (a) the presence of excess zeros using a score test (van den Broek, 1995) (b) the type of excess zeros using a data generation process (Dwivedi et al., 2014b), (c) the presence of overdispersion after removing excess zeros using an overdispersion test (Cameron & Trivedi, 1990), and (d) selection of the best-fit model using model fit criteria (Fernandez & Vatcheva, 2022). Although we followed Fig. 1 for comparing the performance of different models, this step-by-step approach should not be used for prespecifying the primary models in the statistical analysis, but may be used for validating the robustness of findings from different models in sensitivity analyses. The absence of polyps explained by the LR model was considered as no detected polyps, while no detected polyps explained by the count models are considered as undetected polyps. The distribution of polyps was described using the frequency and proportion as well as the average number of polyps per patient with standard deviation (SD) between groups for each dataset. We applied the standard LR, PR, RP, NB, ZIP, ZIRP, ZINB, ZIGP, ZHP, and ZHNB models in all three datasets to compare the presence of detected polyps and the number of detected polyps between (a) CAC and SC procedures in the example dataset 1, (b) segmental and non-segmental protocols in the example dataset 2, and (c) presence and absence of a fellowship program in dataset 3. We compared PDR using the same estimand, by computing marginal effect sizes, including the marginal difference in PDR and the marginal odds ratio of PDR (Long & Jeremy, 2014). Since a bias-corrected Voung test (Desmarais & Harden, 2013) is not appropriate for testing the presence of excess zeros, we primarily used a score test (van den Broek, 1995) followed by a bias-corrected Voung test to determine the presence of excess zeros. A significant p-value of the score test and bias-corrected Vuong test indicates the presence of an excess amount of no polyps. We applied an overdispersion test (Cameron & Trivedi, 1990) and an overdispersion parameter test in the NB regression (Molla, Muniswamy & Hines, 2012) after removing excess zeros predicted by inflated models to evaluate if there was the presence of an over- or equi-, or underdispersion issue beyond excess zeros. In the presence of excess zeros in all three datasets, the final inflated model was selected based on model fit indices, including the goodness of fit of the model and information criteria. We summarized the performance of all count models using model fit indices such as the Akaike Information Criterion (AIC), Bayesian Information Criterion (BIC), and negative log likelihood (LL). Lower values of AIC and BIC indicate a better fit of the model, while higher negative LL values indicate a better fit among comparative models. We used the chi-square goodness of fit test to assess differences in observed counts vs. expected counts. A non-significant p-value of the goodness of fit test indicates a good fit of the model. Among all fit indices, BIC was primarily used to examine the relative fit of the model across all alternative models. Although we used multiple criteria to provide reasons for using a specific model and the most appropriate model describing polyp data in three illustrative studies, the primary models for analyzing data must be specified a priori through the study registration or statistical analysis plan (SAP) registration, particularly in interventional studies based on theoretical and simulation considerations. The series of tests and model fit indices used in the study are useful for evaluating assumptions required in the application of primary models and conducting additional sensitivity analyses to validate the findings obtained from primary analyses. The purpose of reporting multiple model fit indices and tests to evaluate assumptions was to validate the robustness of the findings, not to create selection bias issues.

In the primary analysis of randomized studies, the models were adjusted for age, sex, and total procedure time and validated by unadjusted analysis (without any covariate adjustment). In the analysis of the observational study, a priori, all the critical covariates, including age, sex, BMI, time of the procedure, and sedation method, were adjusted in the models. Some adjusted variables were removed from the analyses if models, especially ZINB or ZIGP models, did not converge or produce reasonable coefficients. Owing to limited counts of polyps beyond six, we merged six or higher polyp counts in the primary analysis and validated them by considering the actual polyp counts in each dataset. The results of the example dataset 3 were also validated after imputing missing data using the multivariate normal method (Sullivan et al., 2017). In primary analyses, we followed statistical analysis guidelines such as evaluation of effect modification and non-linearity before finalizing models as recommended (Dwivedi, 2022).

The usefulness of inflated or hurdle models was also determined using simulation studies. The purpose of these simulations was to demonstrate how each statistical method might perform using realistic data but differing from the example datasets used in the study. We conducted three simulation studies (a) using the parameters meeting the criteria of the example dataset, (b) considering the least proportion of zeros (no detected polyps) and no presence of sampling zeros, and (c) considering misspecification of the distributional assumption for each dataset example. For each simulation condition, we generated two random variables, including one binary variable following the binomial distribution with a specified probability and the other variable following a count distribution following either Poisson or NB distribution. A product of these two variables yields a zero-inflated count distribution. In the simulation studies, for example, datasets 1 and 2, we did not include a confounder considering a randomized design, while we introduced a continuous normally distributed confounder for example dataset 3, considering the observational study design. We followed simulation experiments using our previous studies (Dwivedi et al., 2014a, 2014b) and computed average bias (difference between estimated coefficient and true coefficient), relative bias (difference between estimated coefficient and true coefficient divided by the true coefficient), 95% confidence width measuring precision in the estimates, and coverage probability of including the true value in the 95% confidence interval on 1,000 replications. The results of LR were summarized with OR, and the results of count models were summarized with RR along with corresponding 95% confidence interval (CI) and two-sided p-value. We also showed the factors associated with detected polyps and the amount of detected polyps obtained using ZIRP and ZINB analyses compared to LR and PR models for the example dataset 3. The estimated PDR from each model was summarized for the primary models. A p-value less than 0.05 was considered statistically significant, and the Stata v17.0 (StataCorp, Lakeway Drive, College Station, TX, USA) package was used for all the statistical analyses.

Results

Discussion

In our study, we observed a significant effect of the CAC in detecting a higher number of polyps compared to SC, without any difference in PDR between procedures using a ZIP model. However, the standard approaches produced potentially inaccurate findings that CAC was neither associated with higher PDR nor with the number of detected polyps. In another example dataset, both the standard LR and ZIP models produced the same findings that there was slightly higher PDR in the segmental withdrawal time protocol compared to the standard withdrawal time protocol. The application of the ZIP model also suggested an insignificant difference in detected polyp counts between the segmental protocol and the non-segmental protocol. In another data example, our proposed analytic approach using ZIRP or ZINB suggested that the number of detected polyps increased after the introduction of the GI fellowship program. The effect sizes (OR/RR) obtained from the inflated models were different than standard approaches, indicating the standard approaches may produce biased estimates, as confirmed by our simulation studies, even if the conclusions do not change. The standard count models may produce biased estimates in the presence of excess zeros, which is the case in our example datasets. The logistic regression may produce an inefficient model (because of data dichotomization), but it may also produce biased estimates, especially in the presence of sampling zeros, which are likely to be observed in colonoscopic studies. While the LR suggested most cofactors were associated with detected polyps, the application of the ZIP model suggested that fellowship was the key factor associated with higher PDR. Contrary to the standard count models, the ZINB model produced factors such as male sex, older age, time of the procedure, and sedation type associated with the polyp counts, while only BMI was associated with PDR. These findings clearly support the application of zero-inflated models over the standard approaches. The application of zero-inflated models provides complete information, including the differences in the presence or absence of polyps as well as differences in the extent of detected polyps compared to single distribution models. Moreover, the marginal estimates obtained from zero-inflated models, particularly ZIP models, yielded more precise estimates for all three datasets, indicating that there is a gain in the use of zero-inflated models in the presence of excess zeros. Many colonoscopy screening trials might have reported potentially inadequate conclusions or inappropriate effect sizes due to the inappropriate use of statistical analysis methods. It is known that inappropriate use of methods, not appropriate with evidence-based biostatistics, may lead to false conclusions (Dwivedi, 2022; Dwivedi & Shukla, 2020).

Although our study clearly demonstrates that there was an excess amount of no polyps in all three datasets indicating that a zero-inflated count model was indeed required to analyze these studies, the simple PR model was sometimes applied to analyze colorectal polyps (Chan, Cohen & Spiegel, 2009; Davies et al., 2023; Hillyer et al., 2014; Liu et al., 2020; Othman et al., 2017). In all three datasets, more than half of the patients had no polyps. Multiple studies (Asadzadeh Aghdaei et al., 2017; Lee et al., 2013) showed that PDR ranges between 20% to 59.9%, indicating a disproportionally higher number of no polyps than detected polyps in colonoscopy studies, attesting to the application of zero-inflated models. Although the NB model produced a better or similar fit based on BIC criteria than the inflated models sometimes, the standard count models yield false findings in mixture data distributions (Campbell, 2021), which is also validated in our simulation studies. For example, in the CAC study, we observed a lower proportion of detected polyps, while increased polyp counts were associated with CAC compared to the SC using inflated models. However, the NB model considered all no polyps and undetected polyps as part of a single count distribution, indicating no difference in the number of detected polyps. Therefore, even if the NB model fits well for the dataset, it should be avoided in the case of an excess amount of zeros (Blasco-Moreno et al., 2019; Fernandez & Vatcheva, 2022). Moreover, sometimes the model fit characteristics do not appropriately guide the proper use of the models (Campbell, 2021). However, reporting model fit characteristics is required in some study objectives (Dwivedi, 2022). Our purpose in reporting multiple model fit characteristics was to assess the robustness of findings using a specific model.

Multiple studies support our findings obtained through ZIP models and refute findings obtained from binary models. In a meta-analysis of 11 studies reporting a comparison of CAC with SC, seven studies did not show any significant improvement in PDR with CAC compared to SC (Mir et al., 2017). In contrast, several studies (Kim et al., 2015; Othman et al., 2017) demonstrated that CAC yielded significantly more polyps than SC, supporting our findings obtained through the ZIP or inflated models, but were missed by the standard PR/NB model. Although increasing withdrawal times were shown to be associated with increased PDR (Haghbin, Zakirkhodjaev & Aziz, 2023), the segmental protocol of spending at least 3 min to the right colon increased PDR but did not show a significant benefit over the standard protocol of withdrawal time, as confirmed in our ZIP models for detection of more polyps. Several studies in separate analyses of dichotomized and count forms of colon polyps confirmed our findings obtained through ZIP or ZINB models that increasing fellows’ participation is associated with increased PDR (Qayed et al., 2017) and the detection of more polyps (Elhanafi et al., 2017; Facciorusso et al., 2020; Lee et al., 2011; Qayed et al., 2017).

We observed a significant dispersion issue in the third example dataset, which we could easily address by performing ZIRP analysis. In fact, we suggest using a robust variance estimation approach in ZIP (ZIRP) or standard PR model (RP) regardless of the dispersion issue, as it can estimate the appropriate variance in under- or slightly overdispersion situations (Campbell, 2021; Payne et al., 2018). If there is no overdispersion and no excess zero issues then one should prefer using QP or GP over the ordinary PR model (Harris, Yang & Hardin, 2012). However, if there are no excess zeros but there is an overdispersion then analysts may prefer using the NB model (Fernandez & Vatcheva, 2022). If there is a non-equidispersion issue in addition to excess zeros then it is recommended to use ZINB or ZIGP (Ismail & Jemain, 2007; Yang, Hardin & Addy, 2009). However, these models may produce an overfit issue by yielding unusually high effect sizes if there is an underdispersion in count values, excluding excess zeros (Blasco-Moreno et al., 2019; Campbell, 2021; Fernandez & Vatcheva, 2022) which is the case mostly occurring in colorectal polyp data as we observed in our randomized studies.

We observed unstable coefficients, especially in the analysis of example dataset 2, and observed convergence issues, especially with ZINB and ZIGB for analyzing datasets 1 and 3 without adjusting for covariates. The methodologists need to consider event per variable (EPV) which typically requires 10 EPV for the LR part of the inflated model (Peduzzi et al., 1996), avoid using NB or ZINB models if there is no overdispersion beyond excess zeros (Fernandez & Vatcheva, 2022), avoid adjusting covariates that are not needed for the inflated part of the model, and consider the appropriate sample size (Williamson et al., 2022). In addition, ZHNB may be preferred over ZINB or ZIGP models to avoid convergence issues if there are limited (<10%) sampling zeros. Although we observed similar conclusions with any inflated or hurdle models in both randomized studies, the overdispersion tests showed no significant presence of overdispersion, indicating the use of ZIRP over the ZINB, ZIGP, ZHP, or ZHNB model. We recommend using ZIRP in such a situation, which can estimate robust variance even with a misspecification of the distributional assumption, as shown in different simulation scenarios. However, in the presence of overdispersion along with excess zeros, as we observed in the example dataset 3, we recommend exploring ZINB or ZHNB models along with the ZIRP model. Our simulation studies also confirmed these recommendations. In colonoscopic procedures, there is a chance of having undetected polyps (false or sampling zeros) even if there is an existing polyp (Pamudurthy, Lodhia & Konda, 2020). Therefore, zero-inflated models are the best choice for analyzing colorectal polyps according to the best biostatistics and evidence-based biostatistics practices (Dwivedi, 2022; Dwivedi & Shukla, 2020). If the analysts are unsure about the type of zeros, require ease of interpretation, or face convergence issues or instability in estimated parameters, particularly with ZINB or ZIGP models, then zero hurdle models may be used over the zero-inflated models. In our study, we did not notice differences in hurdle and inflated models for randomized studies. However, considering less precise estimates with hurdle models compared to inflated models as obtained in simulation studies, inflated models, preferably ZIRP should be used for analyzing polyp data. However, hurdle models may be preferred over inflated models to avoid convergence or instability issues in studies yielding overdispersion beyond excess zeros.

We also facilitated bias, precision, and power comparisons for inflated or hurdle models over single distribution models in the presence of excess zeros. However, we did not reevaluate other published trials in this area to demonstrate further change in the conclusions by using inflated models. Another limitation of our study is missing some confounders in dataset 3 and the application of these models in studies with repeated measures design. The major clinical heterogeneity occurs due to the bowel preparation quality and the experience of the gastroenterologists performing the colonoscopy. In randomized studies, these were balanced due to strict randomization. We adjusted these sources of variation in the analysis of the observational study. Although our goal of the study was to demonstrate the usefulness of inflated/hurdle models for analyzing count data in association and interventional studies, there are emerging machine learning models for handling overdispersion and excess zeros which could be compared with standard inflated/hurdle models in prediction studies (Ben Khedher & Yun, 2024; Haghani, Sedehi & Kheiri, 2017; Sidumo, Sonono & Takaidza, 2024). We highlight that colonoscopy studies should use polyp counts as the quality measure, and methodologists should jointly analyze the PDR with polyp counts using inflated models. To the best of our knowledge, our study for the first time demonstrates that count outcome data generated in colposcopy studies should be analyzed with inflated models as opposed to the standard count or binary LR models. We also demonstrate that ZIRP is more appropriate for underdispersed or equidispersed data, which is typically observed in randomized studies. Our work provides a guidance document for selecting an appropriate choice of count models in view of excess zeros, dispersion, and type of excess zeros, and related statistical procedures for analyzing polyp data. For promoting value-based biostatistics practice, we strongly recommend the use of zero-inflated models, particularly ZIRP for detecting the number of polyps in the primary analysis, regardless of checking the fitting parameters.

Conclusions

Studies on colorectal polyps are predominantly analyzed with the LR model, leading to inefficient estimates and potentially inaccurate findings. We observed that colon cancer screening studies yielded an excess amount of no polyps, which were unexplained by the standard count models. We demonstrated that the standard analytic approaches for analyzing colorectal polyps produced biased effect sizes due to poor model fit characteristics, inaccurate interpretations, and conclusions as compared to the inflated models. ZIP with robust variance estimation was found to be the optimal method for analyzing polyp data in randomized studies, while the ZINB model showed a better fit in an observational study. Considering the data generation process of colonoscopic studies and simulation studies, we strongly endorse the use of zero-inflated models for evaluating colorectal polyps, regardless of fitting characteristics in colonoscopy screening studies, for proper clinical interpretation of data and accurate reporting of findings.

Supplemental Information