Examining the causal association between 25-hydroxyvitamin D and caries in children and adults: a two-sample Mendelian randomization approach

Background: Prior observational studies have reported that higher levels of vitamin D are associated with decreased caries risk in children. However, these studies are prone to bias and confounding so do not provide causal inference. Genetic variants associated with a risk factor of interest can be used as proxies, in a Mendelian randomization (MR) analysis, to test for causal association with an outcome. The objective was to estimate the causal association between serum 25-hydroxyvitamin D (25(OH)D) (the commonly measured vitamin D metabolite in blood) and dental caries using a two-sample MR approach which estimates the causal effect of an exposure on an outcome. Methods: A total of 79 genetic variants reliably associated with 25(OH)D were identified from genome-wide association studies and used as a proxy measure of 25(OH)D. The association of this proxy measure with three outcome measures was tested; specifically: caries in primary teeth (n=17,035, aged 3-12 years), caries in permanent teeth in childhood and adolescence (n=13,386, aged 6-18 years), and caries severity in adulthood proxied by decayed, missing and filled tooth surfaces (DMFS) counts (n=26,792, aged 18-93 years). Results: The estimated causal effect of a one standard deviation increase in natural log-transformed 25(OH)D could be summarized as an odds ratio of 1.06 (95%CI: 0.81, 1.31; P=0.66) for caries in primary teeth and 1.00 (95%CI: 0.76, 1.23; P=0.97) for caries in permanent teeth in childhood and adolescence. In adults, the estimated casual effect of a one standard deviation increase in natural log-transformed 25(OH)D was 0.31 fewer affected tooth surfaces (95%CI: from 1.81 fewer DMFS to 1.19 more DMFS; P=0.68) Conclusions: The MR-derived effect estimates for these three measures are small in magnitude with wide confidence intervals and do not provide evidence for a causal relationship between 25(OH)D and dental caries.


Introduction
Dental caries is a disease process which can lead to irreversible damage to tooth tissues. Initially hydroxyapatite crystals in the enamel, dentine and cementum tissues are demineralized when acidic by-products from bacterial fermentation of simple carbohydrates lead to low pH and mineral undersaturation in the tooth surrounding fluids. Eventually demineralization is followed by a proteolytic destruction of the organic substances of the tooth tissues and a cavity is formed (Selwitz et al., 2007).
Both genetic and environmental risk factors influence dental caries and there is a need to identify modifiable risk factors which could be targets for effective interventions. Vitamin D has been suggested as a potential modifiable risk factor. There is an inverse association between serum 25-hydroxyvitamin D (25(OH)D) (the commonly measured vitamin D metabolite in blood) and caries in childhood (Kim et al., 2018;Schroth et al., 2016), and potential mechanisms, including tooth mineralization and antibacterial effects, have been suggested. Vitamin D stimulates absorption of calcium (Veldurthy et al., 2016), and phosphate (Fukumoto, 2014), so may be relevant to hydroxyapatite crystal structure and mineralization. Vitamin D induces genomic effects in odontoblasts (dentine formation) and ameloblasts (enamel formation) through vitamin D receptor signalling (Zhang et al., 2008). A 6-year follow up doubleblind randomized clinical trial found that high-dose vitamin D supplementation during pregnancy was associated with reduced odds of enamel defects in the offspring (Nørrisgaard et al., 2019). Therefore, vitamin D deficiencies during tooth formation may cause the hard tissues of the tooth to be more sensitive to demineralization. Vitamin D stimulates the production of antimicrobial peptides, such as cathelicidins, which are effective against opportunistic gram-positive and gram-negative bacteria in the tooth biofilm (Youssef et al., 2011). In addition, insufficient levels of vitamin D have been reported to be associated with atrophy of the salivary glands with impaired saliva secretion and increased caries risk (Scardina & Messina, 2012).
These studies do not allow for direct causal inference but have led to the hypothesis that vitamin D supplementation in childhood may be a way to prevent dental caries at a population level. To test this hypothesis, controlled trials have been conducted. A meta-analysis of controlled clinical trials supported causal effects of vitamin D, but also highlighted substantial heterogeneity in estimates, risk of bias in individual studies and strong evidence for publication bias in the available literature (Hujoel, 2013). There is therefore a need for additional sources of evidence to help clarify the causal effect of vitamin D in dental caries.
Mendelian randomization (MR) is an alternative way to estimate causal effects, when randomized controlled trials are unfeasible or inconclusive. This method uses genetic variation, which is reliably associated with the exposure of interest, as a proxy measure of the exposure. These genetic variations are often single base pair changes in germline genotypes, termed single nucleotide polymorphisms or SNPs. A higher number of 25(OH)D increasing variants is associated with higher average serum 25(OH)D concentration. Due to the essentially random assortment of alleles during meiosis, the proxy measure of 25(OH)D concentration is unrelated to traditional confounding factors at population level. For example, individuals with two 25(OH)D increasing alleles at a particular SNP will smoke no more or less than those who have zero 25(OH)D increasing alleles at the same SNP. This is in contrast to serum 25(OH)D measures which, at a population level, show a strong inverse association with smoking (Kassi et al., 2015).
In addition, since disease processes (caries in this case) cannot alter germline genotypes, the direction of causation is from the proxy measure of exposure to the outcome, therefore reducing the risk of bias from reverse causation (Lawlor et al., 2008). MR has been increasingly used over the past decade, to provide more robust causal effect estimates for a range of risk factors and health outcomes (Davies et al., 2018).
A previous study has assessed the causal effect of 25(OH)D on caries using MR, but was underpowered to provide precise estimates or fully interrogate the assumptions of the MR method This version includes updates made according to the reviewers' feedback, including: • Change in title name to be more informative and include the method of two-sample Mendelian randomization; • Addition of Figure 1, a diagrammatic representation of the stages involved in two-sample Mendelian randomization analysis (methods); • Elaboration on how the variance explained (4.4%) was calculated, using the formula: variance explained ≈2β 2 ƒ(1-ƒ), where β and ƒ denote the effect estimate and the effect allele frequency of the allele on a standardized phenotype, respectively (methods); • Additional sensitivity analyses using 7 independent SNPs that have biological relevance to 25(OH)D to avoid horizontal pleiotropy and confirm the null associations (results); • Chromosome positions have been added to Table 2 and the order of the SNPs has been rearranged in the same way as Table 1 for easier reading (results); • Renaming of the column in Table 3 to "odds ratio/transformed effect" for further clarity; and • Removal of Table 4 and the post hoc power analysis in the methods, results and discussion, as this analysis is redundant given the reported confidence intervals. Instead, the wide confidence intervals have been used to indicate the precision of the observed estimates (discussion).

REVISED
create the opportunity to re-examine the association of 25(OH)D and dental caries.
The objective of this study is to use a two-sample MR analysis to assess the causal role of serum 25(OH)D on three caries related traits (caries in the primary dentition, caries in the permanent dentition in children and adolescents, and caries severity in adults) using data from published genome-wide association studies (GWASs).

Methods
The data for SNP-exposure and SNP-outcomes were extracted from published GWASs. GWASs are used in genetic research to identify genetic variants that are associated with a specific trait. Figure 1 shows the steps involved in performing a two-sample MR.
SNP-exposure: identifying SNPs as 25(OH)D proxies SNPs to be used as proxies for 25(OH)D exposure were identified as those that a) were strongly associated with 25(OH)D exposure (defined as P≤5×10 -8 in at least one published GWAS); b) were identified in a population of European ancestry participants; c) were conditionally independent of other SNPs in the same genomic locus in conditional analysis (Manousaki et al., 2020); and d) had a minor allele frequency of 0.05 (5%) or above in the 25(OH)D GWAS and the caries outcome GWASs. For each variant, information about the strength and direction of association with 25(OH)D was extracted and recorded alongside other information such as variant name and nearest gene.
SNP-outcome: obtaining estimates of the association between 25(OH)D-proxying SNPs and caries outcomes GWASs test for association between millions of SNPs across the genome and diseases such as dental caries. Existing GWAS studies can therefore be used to obtain estimates of the association between genetic proxies for an exposure of interest and an outcome. For this investigation, the estimated effects of 25(OH)D-proxying SNPs on dental caries were obtained by extracting association statistics for these variants from published results.
SNP-caries association estimates for caries in the primary dentition in children and caries in the permanent dentition in adolescents were obtained from a published meta-analysis. This meta-analysis was originally carried out in a consortium, including cohort studies in the USA, UK, Denmark, Finland, the Netherlands, Germany and Australia (Haworth et al., 2018). Each of these studies classified participants as caries-free or caries-affected using clinical examination, index linkage to pre-existing dental records or using intra-oral photographs. One contributing study used child-and parent-reported questionnaires to classify children as caries-free or caries-affected. The relationship between GWAS derived SNPs and caries status was estimated in a logistic regression framework accounting for co-variables such as age, sex and genetic principal components, which aim to control for variation in genetic data which is due to genetic ancestry. The log-odds ratio (OR) effect estimates were combined in a fixed-effects genome-wide meta-analysis. The principal analysis included all studies with phenotypic data, and a sensitivity analysis was performed which excluded participants with questionnaire-derived caries status. Results from the principal meta-analysis were download from the University of Bristol research data repository and are available at https://doi.org/10.5523/bris.pkqcnil6e9ju2nyreblt3mvwf (Haworth & Timpson, 2018). Downladed data were used to extract the SNP-outcome information about the vitamin D-associated SNPs, for caries in the primary and caries in the permanent dentition in children and adolescents, respectively.
SNP-caries association estimates for decayed, missing and filled tooth surfaces (DMFS) were obtained from a genome-wide meta analysis of adult cohort studies based in the USA, Germany, Sweden and Finland (Shungin et al., 2019). In each study, DMFS scores were calculated excluding third molar teeth and root caries, and were derived from surface-level dental charts, which were either obtained as part of the study protocol or obtained via index linkage to pre-existing dental records. Within each study, DMFS scores were residualized on covariates including age, age squared and sex, then residuals were standardized to a mean of 0 and standard deviation of 1. Subsequently, the relationship between SNPs and transformed DMFS scores was estimated using a linear regression framework, and beta coefficients were combined across studies using a genome-wide fixed-effects meta-analysis. The principal analysis included all studies, while a sensitivity analysis was performed excluding one study of Hispanic/Latino participants. Results of the primary and sensitivity meta-analyses were downloaded from the University of Bristol research data repository at https://doi. org/10.5523/bris.2j2rqgzedxlq02oqbb4vmycnc2 (Haworth, 2019). This dataset was used to extract the SNP-outcome information about the 25(OH)D-associated SNPs for the caries severity in adults outcome.

Causal effect estimation
Intuitively, if SNPs which are proxies for 25(OH)D are associated with caries outcomes then there must be an effect of 25(OH)D on caries. This intuition is formalized using a series of statistical tests under an analytical paradigm termed two-sample MR (Burgess et al., 2015). In this investigation, the "TwoSampleMR" package version 0.5.4 in R version 3.4.3 (Hemani et al., 2016) was used. For the primary analysis, causal effect estimates were obtained separately for each SNP and then combined using inverse variance weighted (IVW) meta-analyses to estimate the overall causal effect of serum 25(OH)D on each of the three outcomes.
For caries in the primary dentition and permanent dentition in children and adolescents each causal effect estimate (β xy ) and standard error (SE) from the two-sample MR was converted into an interpretable OR with 95% confidence intervals (CIs) as follows: Leave-one-out sensitivity analysis For each outcome, a leave-one-out sensitivity analysis was undertaken. In this analysis, the IVW analysis was carried out iteratively omitting one 25(OH)D proxying SNP in turn. This analysis aims to identify SNPs with outlying causal effect estimates and ensures that the overall estimate of the causal effect of 25(OH)D on dental caries was not driven by a single or few SNPs with large effects. In cases where a small number of SNPs have large effects on caries (but not on 25(OH)D), this may imply that these SNPs act on caries through mechanisms unrelated to vitamin D (pleiotropic effects) and violate the core assumptions of the MR method.

Weighted median sensitivity analysis
For each outcome, a weighted median (WM) sensitivity analysis of the SNPs was undertaken. The WM provides a consistent estimate of the causal effect, provided that 50% of the weight in the analysis stems from non-pleiotropic, and therefore valid, variants.

MR-Egger regression sensitivity analysis
For each outcome, a MR-Egger regression of the SNP-outcome on the SNP-exposure with the y-intercept unconstrained was carried out. The y-intercept tests for the presence of directional pleiotropy, since when the SNP-exposure association is zero the SNP-outcome association should also be zero. The MR-Egger regressions rely on the instrument strength independent of direct effect (InSIDE) assumption that the strength of the SNP-exposure association should not correlate with the strength of any pleiotropic effects across the instrumental variants (Bowden et al., 2015).

Strength of instrumental variables
Instrumental variable (IV) estimates can suffer from weak instrument bias, which arises when confounders in the genotypic subgroups in samples are not perfectly balanced. If the IVs are weak, they explain less variation in the phenotype and therefore the difference in confounders between the subgroups might explain more of the variation in phenotype. If this bias is present, a confounded false-positive association between exposure and outcome may be found (Burgess et al., 2011). The total variance explained in 25(OH)D by the SNPs was estimated as the sum of the variance explained by each SNP in the experiment. The variance of each SNP was given using the formula provided by the authors of the 25(OH)D GWAS paper: variance explained ≈ 2β 2 ƒ(1-ƒ), where β and ƒ denote the effect estimate and the effect allele frequency of the allele on a standardized phenotype, respectively (Manousaki et al., 2020).

Ethical approval
This analysis of previously published data was conducted in accordance with principles described in the Helsinki declaration and all addition requirements within the United Kingdom. All participants participating in the published studies which contributed to this analysis gave informed consent, as described in the respective publications.

SNP-exposure: identifying SNPs as 25(OH)D proxies
In total, 83 SNPs met the criteria for inclusion as proxies for 25(OH)D; these were selected from the most recent published GWAS of 25(OH)D (Manousaki et al., 2020). However, seven SNPs associated with 25(OH)D identified from the SNPexposure stage were not present in the outcome data for all three traits and one SNP (rs200454003) was not present in the outcome data for caries in primary teeth and caries in permanent teeth in paediatric populations. Proxy SNPs with the highest linkage disequilibrium in European populations (r 2 ) to the seven missing SNPs were identified. Proxy SNPs were used for four SNPs (rs2934744, rs7650253, rs3822868, rs201501563) where the r 2 was 0.7 or greater. The remaining three missing SNPs (rs145432346, rs200641845, rs3775150) were excluded from the exposure data as no suitable proxy was identified, leaving 79 SNPs for analysis for caries in primary teeth and caries in permanent teeth in paediatric populations.
For DMFS, two SNPs (rs10127775 and rs10832289) were removed for being palindromic with intermediate allele frequencies (to harmonize the data so that the effect of the variants both exposure and outcome corresponded to the same allele), leaving 78 SNPs for analysis. Information about the 79 SNPs is provided in Table 1.

SNP-outcome: obtaining estimates of the association between 25(OH)D-proxying SNPs and caries outcomes
For caries in primary teeth, binary data for 17,035 children (aged 3-12 years) were available from nine studies of European  (Table 2).

Causal effect estimation
In paediatric populations, the estimated causal effect of 25(OH)D on caries in primary and permanent teeth was essentially null (Table 3). For caries in primary teeth the estimated effect was small in magnitude with CIs crossing the null Table 2. Summary of SNP-outcome association statistics for the three caries outcomes: caries in primary teeth, caries in permanent teeth and decayed, missing, and filled surfaces (DMFS). In adult populations, the estimated causal effect of 25(OH)D was in the direction of a small protective effect, with higher 25(OH)D associated with fewer caries-affected tooth surfaces (0.31 fewer caries-affected tooth surfaces [95%CI: from 1.81 fewer DMFS to 1.19 more DMFS, P= 0.68). This effect was again small in magnitude and the CIs did not provide evidence supporting a causal relationship (Table 3).

Leave-one-out sensitivity analysis
Leave-one-out sensitivity analysis suggested that no single SNP was strongly driving the IVW point estimates for each trait, since there were no cases where excluding one SNP resulted in dramatic changes in the overall result. In addition, all the CIs overlapped with the other causal estimates in each sensitivity analysis, further suggesting that all SNPs estimated a single common exposure (Figure 2-Figure 4).

Weighted median sensitivity analysis
Weighted Median sensitivity analysis for all three traits provided effect estimates similar in magnitude to the IVW primary analysis (Table 3). This suggests that the 25(OH)D variants were valid instruments.

MR-Egger sensitivity analysis
The MR-Egger sensitivity analysis for all three traits provided effect estimates similar in magnitude to the IVW primary analysis and WM analysis ( Sensitivity analysis using 7 independent SNPS Sensitivity analyses using 7 independent SNPs, located in the the GC and CYP2R1 loci, that have biological relevance to 25(OH)D was carried out, to avoid horizontal pleiotropy and confirm the null associations. Results of this sensitivity analysis were consistent with the principal analysis (caries in primary

Strength of instrumental variables
The 81 SNPs combined were considered to be a strong proxy for 25(OH)D in the MR study. The SNPs explain approximately 4.4% of the total variation in 25(OH)D levels in populations of European ancestry.

Discussion
This study tested for causal effects of serum 25(OH)D on dental caries using an MR approach. The main findings are near-null causal effect estimates for all three dental caries traits which do not provide evidence to support a causal association between 25(OH)D and dental caries. The effect sizes seen in this study agree with some reported associations in the observational literature, for example a weak association with a similar OR (0.97) found by Gyll et al. (2018) but disagrees with other One limitation of the MR approach is statistical power, since this method requires large samples to produce precice causal estimates. The apparent disagreement in the results of this study and previous literature might be due to statistical power, the wide 95% CIs for all three caries outcome traits reflect the need for much larger sample sizes. We did not see large effect sizes which could be interpreted as suggesting that confounding in the observational literature and publication bias in results of controlled clinical trials has potentially over-estimated the effect sizes, however we also consider other explanations as discussed below.
To provide valid inference, the MR method makes three assumptions. First, the genetic variants used as proxies for the exposure need to be robustly associated with the exposure of interest. This assumption has been satisfied in the current experiment since the variants were identified from published GWAS as strongly associated with 25(OH)D (P values all ≤ 5×10 -8 ).
In addition, the F-statistics are >10, therefore the IVs are considered to be strong, owing to the large sample size of available outcomes for the three traits under study (Burgess et al., 2011). It therefore appears unlikely that weak instrument bias explains the lack of causal association seen in the present study.
The second assumption states that the genetic variants must be independent of confounders of the exposure-outcome association. In general, this assumption is difficult to test in two-sample MR experiments, but violations might arise due to population stratification (different populations inheriting haplotype blocks in different frequencies due to differences in ancestries). To help address this, nearly all participants included in the MR analysis were of European ancestry and individuals are assumed to only differ with respect to the 25(OH)D loci under study. In addition, the SNP-25(OH)D SNP-caries association statistics were obtained from models which adjusted for population substructure, as described in the respective GWAS publications. There was little evidence for inflationary bias in the summary statistics of these models, suggesting that population stratification was well-controlled. Finally, sensitivity analyses were performed using more stringent exclusion criteria for ancestry in adults, and the results of these analyses agreed with the primary analysis. These approaches help protect against violations of the second assumption, but as a result of sample restriction, the results from this study are only applicable to individuals of European ancestry and are not necessarily generalizable to other ethnic populations. As a more general observation on external validity, the children included in this study were recruited from wealthy countries with longstanding public health messages regarding both vitamin D and dental caries. Thus, the prevalence of both vitamin D deficiencies and caries may be lower in the study population than in other groups, and care is needed in extrapolating the findings to populations in other countries or indeed historical studies in the same countries.
The third assumption states that the genetic instrument must be associated with dental caries, only through its effects on 25(OH)D and not via any alternative pathways (referred to as pleiotropy). If any single SNP had strongly pleiotropic effects acting on the outcome through a mechanism other than through 25(OH)D, this would be shown in the leave-one-out sensitivity analysis as an outlier. No such outliers were identified. Groups of SNPs with pleiotropic effects would result in a discrepancy between the results of IVW and WM analysis, which we did not observe. Finally, groups of SNPs with unbalanced horizontal pleiotropic effects would result in a detectable intercept term in MR-Egger regression analysis.
This experiment used a two-sample design which places natural limits on the scope of the experiment. For example, it was not possible to test for non-linear effects of 25(OH)D or describe the full characteristics of the SNP instruments. However, within these limits we have tested for we did not find evidence of violations of the assumptions of the method.
In summary, MR has provided a genetic approach to assess causality between 25(OH)D and dental caries, an association which is currently not fully understood. This study did not find evidence supporting a clinically relevant causal association. Although the assumptions required by this method appeared to be valid, results were similar across different caries outcomes and findings were consistent in sensitivity analysis, we acknowledge that statistical power was a limitation. At this moment in time the results do not suggest that vitamin D supplementation is likely to be an effective population-level risk reduction strategy for alleviating the burden of disease from dental caries where there is no suspicion of vitamin D insufficiency. In the future, larger GWAS for caries may provide more precise quantification of the role of vitamin D or other modifiable risk factors in the aetiology of this complex and important disease.

Robert J Schroth
Rady Faculty of Health Sciences, University of Manitoba, Winnipeg, Canada Thank you very much for the opportunity to review this manuscript. This is a very interesting and novel study exploring whether there is evidence of a causal relationship between vitamin D and caries in children and adults. The authors have used a Mendelian Randomization (MR) approach to estimate the causal effect between serum 25(OH)D and caries.
Title: This reviewer would suggest a change in the title so that it would be more informative to readers and indicate the type of study it is. Perhaps it would be better as "Examining the causal relationship between 25-hydroxyvitamin D and caries in children and adults: a MR analysis approach." Abstract: the use of the term "null hypothesis" may not be familiar for readers from some regions of the world. It would be more clear to just state what the hypothesis was and whether the findings support the causal relationship or not.
Introduction: is appropriate. The only comment would be that the authors have reported on findings from cross-sectional studies. Is it worth mentioning findings from any prospective studies looking at 25(OH)D during periods of tooth formation and caries at a later time point?
Methods: perhaps the addition of an extra sentence or two about genome-wide association studies (GWASs) would help future readers. Would it be possible to add a small figure that shows the various GWASs whether data were obtained from to conduct this MR study?
Results: well described.
Discussion: as recommended earlier, it might be best to rephrase and avoid the use of the term "null hypothesis".
Tables and figures: appropriate.

Is the work clearly and accurately presented and does it cite the current literature?
Thank you for the suggestion, the term null hypothesis has been removed and reworded. Thank you for the comment. We have added information from an additional study carried out during tooth formation to the introduction. Thank you for the comment. A new figure (Figure 1) has been added to the methods section to further explain the two-sample MR method.
"The data for SNP-exposure and SNP-outcomes were extracted from published GWASs. GWASs are used in genetic research to identify genetic variants that are associated with a specific trait. Figure 1 shows the steps involved in performing a two-sample MR." " Figure 1 5. Discussion: as recommended earlier, it might be best to rephrase and avoid the use of the term "null hypothesis".
Thank you for the comment the term null hypothesis has been removed from the discussion as well.
"The main findings are near-null causal effect estimates for all three dental caries traits which do not provide evidence to support a causal association between 25(OH)D and dental caries." I would like to discuss the necessity for post hoc power analysis with the authors. The power analysis that the authors applied is designed mainly for one-sample MR 4 . As confidence intervals (CIs) indicate the precision of the observed association estimates and reflect the sample size, post hoc power analysis is generally not recommended 5 .

3.
Minor comments: Please write clearly that the study is a two-sample MR in the title, abstract and the introduction section. 1.

In the Methods:
One of the selection criteria for SNPs used as proxies for 25(OH)D exposure is "SNPs had independent effects". The used linkage disequilibrium (LD) R 2 cut-off value should be given. 1.
GWAS of DMFS: DMFS as a count variable has a strongly positively skewed distribution with a large stack of zero counts for those without caries. How the zeros were handled when DMFS scores were generated? 2.

In the Results
The F statistic is usually reported in one-sample MR studies to represent strength of the instrumental variables (IVs). In a two-sample setting, the R 2 statistic (a measure of the variance in the exposure explained by the genetic variants) is more proper. In the Results, it reads "In total, these SNPs explain approximately 4.4% of the total variation in 25(OH)D levels in populations of European ancestry". How was this 4.4% calculated? In the study by Manousaki et al., 138 SNPs explained 4.9% of the variation of 25(OH)D level.

1.
For Table 2, could the authors rearrange the order of SNPs in the same way as in Table 1 and give chromosome positions for easier reading? 2.
The tables and figures should be interpretable by themselves. Add more details in table footnotes and figure legends: e.g. Table 1. Number of participants should be given for the SNP-25(OH)D association. According to the study by Manousaki et al., the effect estimates and P values were based on meta-analysis (n=443,734).
1. Table 3: use OR instead of "Transformed effect", OR and beta correspond to genetically determined one standard deviation increase in natural log-transformed 25(OH)D.

3.
The existence of pleiotropy likely changes the estimates away from the null association. However, in theory it is possible that pleiotropy leads to a false negative finding. It would be advisable if the authors could perform extra sensitivity analyses using only independent SNPs related to four loci (GC, NADSYN1/DHCR7, CYP2R1, CYP24A1 ) that have biological relevance to 25(OH)D to avoid horizontal pleiotropy and confirm the null associations.
Thank you for this comment. As suggested we have included additional sensitivity analysis using a subset of variants (rs705117, rs222026, rs201501563, rs186881826, rs11723621, rs10832289 and rs10832218). These are located in the GC and CYP2R1 loci. These variants were selected as they were available in the outcome data and they had a minor allele frequency of 0.05 or above as stated in the SNP-exposure methods. The results of this analysis were concordant with the main analysis and had similar interpretation. We have added a paragraph to the results section: "Sensitivity analysis using 7 independent SNPS Sensitivity analyses using 7 independent SNPs, located in the the GC  2. I would like to discuss the necessity for post hoc power analysis with the authors. The power analysis that the authors applied is designed mainly for one-sample MR. As confidence intervals (CIs) indicate the precision of the observed association estimates and reflect the sample size, post hoc power analysis is generally not recommended.
Thank you for this comment. We agree that the power analysis applied is redundant given the reported confidence intervals and have therefore removed Table 4 and this analysis from the methods, results and discussion. The discussion has been reworded to describe the wide confidence intervals as an indication of the greater sample size that is needed to provide more precise effect estimates.
"The apparent disagreement in the results of this study and previous literature might be due to statistical power, the wide 95% CIs for all three caries outcome traits reflect the need for much larger sample sizes. We did not see large effect sizes which could be interpreted as suggesting that confounding in the observational literature and publication bias in results of controlled clinical trials has potentially over-estimated the effect sizes, however we also consider other explanations as discussed below."

Please write clearly that the study is a two-sample MR in the title, abstract and the introduction section.
Thank you for this comment. This description has been added to the title, abstract and introduction section. 4. One of the selection criteria for SNPs used as proxies for 25(OH)D exposure is "SNPs had independent effects". The used linkage disequilibrium (LD) R 2 cut-off value should be given.
Thank you for this comment. The SNPs were chosen based on a conditional analysis implemented in GCTA, as reported in the original manuscript. Variants were selected if they reached genome-wide significance conditional on the lead variant in the same locus. Variants situated more than 20,000 base pairs away and variants with LD R 2 >0.9 were excluded. For further details of GCTA-COJO analysis please refer to the section "Approximate Conditional Association Analysis" in the following paper titled "Genome-wide Association Study for Vitamin D Levels Reveals 69 Independent Loci": https://pubmed.ncbi.nlm.nih.gov/32059762/.
"c) were conditionally independent of other SNPs in the same genomic locus in conditional analysis (Manousaki et al., 2020)." 5. GWAS of DMFS: DMFS as a count variable has a strongly positively skewed distribution with a large stack of zero counts for those without caries. How the zeros were handled when DMFS scores were generated?
Thank you for this comment. The GWAS for DMFS was carried out in adult populations. DMFS scores were regressed on age, age squared, and study specific-covariates and the residuals were plotted. We did observe skewed distributions for other traits (for example number of teeth), but not for DMFS, where the residuals were approximately normally distributed. This may be because the study population included middle aged and older adults with high disease burden. DMFS residuals were subsequently treated as a continuous trait in GWAS analysis.
6. The F statistic is usually reported in one-sample MR studies to represent strength of the instrumental variables (IVs). In a two-sample setting, the R 2 statistic (a measure of the variance in the exposure explained by the genetic variants) is more proper. In the Results, it reads "In total, these SNPs explain approximately 4.4% of the total variation in 25(OH)D levels in populations of European ancestry". How was this 4.4% calculated? In the study by Manousaki et al., 138 SNPs explained 4.9% of the variation of 25(OH)D level.
Thank you for this comment. While we aimed to provide some information on the potential strength of the MR instrument, we accept that the F statistic is potentially confusing and have therefore removed this part. We have amended the Results section: The R 2 of the 25(OH)D reported in the results was estimated as the sum of variance explained by each SNP in the experiment. The variance of each SNP was given using the formula provided by the authors of the Vitamin D GWAS paper; variance explained ≈2β 2 ƒ(1-ƒ), where β and ƒ denote the effect estimate and the effect allele frequency of the allele on a standardized phenotype, respectively. This information has been added to the methods section.
"The total variance explained in 25(OH)D by the SNPs was estimated as the sum of the variance explained by each SNP in the experiment. The variance of each SNP was given using the formula provided by the authors of the 25(OH)D GWAS paper: variance explained ≈ 2β 2 ƒ (1-ƒ), where β and ƒ denote the effect estimate and the effect allele frequency of the allele on a standardized phenotype, respectively (Manousaki et al., 2020)." Table 2, could the authors rearrange the order of SNPs in the same way as in Table 1 and give chromosome positions for easier reading?

For
Thank you for the suggestion, Table 2 has been rearranged to the same way as Table 1. Thank you for the comment. More details have been added to the footnotes and figure legends. Table 3: use OR instead of "Transformed effect", OR and beta correspond to genetically determined one standard deviation increase in natural log-transformed 25(OH)D.

9.
Thank you for the comment. The DMFS effect estimate is not on a odds ratio scale and therefore we prefer to use the term 'transformed effect' for this variable. The column has been renamed to "Odds ratio/transformed effect". Thank you for this comment. We have added the information as suggested.
Competing Interests: No competing interests were disclosed.