Broad Spectrum Epidemiological Contribution of Cannabis, Tobacco and Alcohol to the Teratological Prole of Northern New South Wales: Geospatial and Causal Inference Analysis

Background. Whilst cannabis commercialization is occurring rapidly guided by highly individualistic public narratives, evidence that all congenital anomalies (CA) increase alongside cannabis use in Canada, a link with 21 CA’s in Hawaii, and rising CA’s in Colorado indicate that transgenerational effects can be signicant and impact public health. It was therefore important to study Northern New South Wales (NNSW) where cannabis use is high. Methods. Design: Cohort. 2008–2015. NNSW and Queensland (QLD), Australia. Participants. Whole populations. Exposures. Tobacco, alcohol, cannabis. National Household Surveys 2010, 2013. Main Outcomes. CA NNSW-QLD comparisons. Geospatial and causal regression. intestinal atresia), body wall (gastroschisis, diaphragmatic hernia) and other (hypospadias) (AVTPCDSGDH) CA’s. In linear modelling cannabis use was signicantly linked with anal stenosis, congenital hydrocephalus and Turner syndrome (ACT) and was signicantly linked in borderline signicant models (model P < 0.1) with microtia, microphthalmia, and transposition of the great vessels. At robust and mixed effects inverse probability weighted multivariable regression cannabis was related to 18 defects. E-Values in spatial models were generally > 1.3 ranging up to 3.8 × 10 30 making uncontrolled confounding unlikely.

A recent report on Canada demonstrated that total congenital defects were three times more common in the northern Territories which smoked more cannabis than the Provinces and that the association was robust to socioeconomic adjustment [2]. A recent study from Colorado across the period of cannabis legalization showed that many defects rose parallel to increased cannabis consumption including all chromosomal defects (ACD), Downs syndrome and several cardiovascular defects including atrial septal defect (ASD) and patent ductus arteriosus (PDA), common defects which had not been previously linked with prenatal cannabis exposure (PCE) [3]. It was calculated that in Colorado over 11,000 extra defects occurred 2000-2014 related to increased cannabis use [3]. An Hawaiian study found that 21 defects were increased in mothers who were exposed only to cannabis [4].
Whilst some of these studies have used sophisticated geospatial modelling techniques [2] all epidemiological research is fundamentally associational in nature. However similar ndings elsewhere strengthens the evidence base.
Northern New South Wales (NNSW, NSW) is a well known drug using and cannabis cultivation area 760km from Sydney but only 180km from tertiary pediatric care centres in Brisbane and 111km from Southport both in Queensland (QLD). Although lying within New South Wales administratively many of its neonatal CA's are evacuated to tertiary pediatric hospitals in Queensland under the Neonatal Retrieval Scheme (NRS) [5] and their data thus appears in Queensland statistics. This therefore presents an ideal opportunity to directly compare NNSW and Queensland neonatal epidemiology.
Our hypothesis was that cannabis use would be associated with increased congenital anomalies and was formulated prior to study commencement.

Methods
Data. Data on congenital anomaly rates for Queensland Health service areas including northern New South Wales was taken from the Congenital Anomaly Linked File (CALF) from Queensland Health [6].
Annual data by area has not been publicly released. Data on maternal age was from the QLD and NSW annual Mothers and Babies reports [7,8]. CALF data includes numbers, rates and con dence intervals for the data. Drug use data for last month cigarette use, last month binge alcohol and last year cannabis use by area was obtained from the Australian Institute of Health and Welfare from the National Drug Strategy Household Survey (NDSHS) 2010 and 2013 [9] and averaged to obtain a mean rate by area across this period pursuant to our custom data request. Data was matched manually between drug use and congenital anomaly datasets. Areal shape les were taken from the Australian Government national website [10]. The northern coastal area of NSW was added on to the Queensland Health shape le. This depiction of the NNSW catchment area is illustrative only and not intended to be exact as the geographic boundaries of the NRS are not de ned [5].
Congenital anomalies were de ned as cannabis related or not based on a literature review and recent reports [2, 3, 11-13] particularly [4].
Patient and Public Involvement. Patients were involved in this research from the very outset. Many patients are concerned about reproductive health outcomes. In this area there is signi cant local concern about the implications of widespread drug use in nearby areas. Patients are therefore frequently concerned about optimal long term outcome for their offspring. Patients who are not pregnant or who are not considering starting a family are concerned abut the possible teratogenic outcomes for people in the reproductive age group and the likely impact this might have on the wider community. Research questions in the present study therefore including the outcome measures considered were developed and informed by their priorities experiences, preferences and concerns. Patients and the public are of the view that all available existing datasets relating to this issue should be carefully investigated from this perspective. Patients are concerned also with the outcomes of this study. Our patient are happy to assist with the dissemination of results by means of word of mouth and social media technologies and keen to see such information disseminated widely across the community both locally and internationally.
Statistics. Data was processed in RStudio version 1.2.1335 based on R version 3.6.1 on 16 th April 2020. Two-by-two tables were analyzed in package epiR using epi.2by2. Graphs were drawn with R-Base and ggplot2 and in Excel. Maps were drawn using sf ("simple features") package. Principal component analysis was conducted using the psych package. Linear regression was performed in Base-R. Batch extraction of all linear model coe cients by different defects was performed with broom and purrr packages. Links between neighbouring areas sharing an edge or corner ("queen"-relationships by analogy with chess moves) were derived with the poly2nb function from spdep and edited as indicated.
This neighbourhood map was used to calculate the geospatial weights matrix for spatial regression.
Geospatial regression was performed with the spreml function from package splm [13,14] using the derived spatial weights matrix. All spatial models used a full error structure of Kapoor Kelejian and Prucha [15] and had serially correlated remainder errors and random effects (sem2srre). The appropriateness of this error structure was formally tested by substituting various alternative forms and comparing results including the logLik and spatial Hausman tests. Models were spatially lagged and not lagged as indicated.
Inverse probability weights (IPW) were derived with the ipw package in R using cannabis use as the exposure of interest, tobacco in the numerator and tobacco and alcohol in the denominator. IPW weights were then used in robust regression models conducted in the R package survey, and in mixed effects models in the R package nlme to generate datasets pseudo-randomized for cannabis exposure. This allowed causal relationships to be assessed. E-Values were calculated with the R package EValue to quantify the degree of association some unmeasured confounder would require with both dependent and independent variables to explain away the observed effect.
For all regression models model reduction was by the classical method with sequential deletion of the least signi cant term. Missing data was casewise deleted at multivariable regression. P<0.05 was considered signi cant.

Results
Input data is shown in an online supplementary csv le. Supplementary Table 1 provides comparative congenital anomaly data between QLD and NSW by both numbers and rates including defect relationship to cannabis [6]. Denominator data was calculated from the numbers and rates supplied in that le. It was veri ed from the annual QLD Health Mothers and Babies reports 2008-2015 which show 509,095 births in this period [8]. The "Interstate and Overseas" designation in the CALF le includes offshore islands such as Christmas, Norfolk, Cocos and Lord Howe Islands which together have a population of 4,518. The prime catchment area of the NRS is Northern NSW which has a population of 296,531 [10]. Hence only 1.5% of the population in this designation is likely to come from outside NNSW. The view that the "Interstate and Overseas" designation refers primarily to NNSW is con rmed by QLD Health Ministerial correspondence (Minister Steven Miles, 05/04/2018). The denominator gure calculated for NNSW in this manner is 4,800 births.
It should be noted that NNSW birth defect data also appears in NSW Health records [7]. One notes that the rates of congenital anomalies reported for this region in the NSW Mothers and Babies reports are about half those of the rest of the state. This is presumably related to the relocation of many cases into Queensland through the NRS. Queensland congenital anomaly rates are much higher than those reported elsewhere so it is not possible simply to combine NSW and QLD Health reports. Therefore this report is limited to consideration of the QLD Health CALF le only.
Drug use data is shown in Supplementary Table 2. It is noteworthy that the Richmond-Tweed NNSW area has a middle ranking for tobacco and alcohol use, but a rst ranking for cannabis use.
Maternal age is a major factor bearing on congenital anomaly rates and it is known to be strongly linked with chromosomal anomaly rates. Interestingly CALF Table 1 shows rises in the rates of several defects ( Supplementary Fig. 1) including CVS defects, atrial septal defect (ASD) and ventricular septal defect (VSD) which are highly signi cant (Supplementary Table 3). Intriguingly the mean incidence of daily smoking tobacco and high risk alcohol use dropped across this period and annual cannabis use rose from 10.5-11.3%. The principal component of the combination of cardiovascular, gastrointestinal and respiratory anomalies also rose Page 6/30 signi cantly across this period. These data suggest that cannabis may be a more potent and more important teratogen than tobacco and alcohol.  Figure 1 shows a qualitative choropleth map-graph for the major CA classes. The yellow zones re ect high incidence and dark blue low incidence.
Supplementary Figs. 2-4 present choropleth maps of CA incidence by area. Supplementary Fig. 5 shows chromosomal anomaly incidence for which data is available. Figure 2 was drawn in Excel and shows the con dence intervals from CALF for common, intermediate frequency and rare defects for cannabis-related (CRD) and cannabis not related (CNRD) defects. For most of the cannabis-unrelated defects the con dence intervals overlap. For most of the cannabis-related defects the con dence intervals either do not overlap, or are near the lower end of the QLD C.I.'s. Supplementary Fig. 6 expands this list for rare defects and continues this trend.
Supplementary Fig. 7 compares the QLD and NNSW CA rates.
Supplementary Fig. 8 compares all the rate ratios of defects using the quoted rates in the CALF le. Supplementary Fig. 9 makes a similar comparison with log rates and shows that most of the cannabisrelated defects are more common in NNSW.
CNRD were more common in QLD ( Table 4 lists the PR's, attributable fraction in the exposed (AFE) and attributable fraction in the population (AFP) along with their C.I.'s and applicable P-values for all defects and defect classes. Supplementary Fig. 10 illustrates the PR's and C.I.'s for CRD and CNRD. Figure 3 shows the AFE's and C.I.'s for CRD and CNRD. Supplementary Fig. 11 shows the AFP's and C.I.'s for CRD and CNRD. P-values are illustrated in Supplementary Fig. 12.
Supplementary Fig. 13 shows ve main defect classes charted against the use of tobacco, binge alcohol and cannabis. Rising trends with cannabis seem to apply to CNS, cardiovascular and chromosomal anomalies.
Supplementary Fig. 14 charts all 55 anomalies and anomaly classes against tobacco use.
Supplementary Fig. 15 performs a similar function for binge alcohol.
When a similar exercise is undertaken for cannabis exposure rising trends appear in several defects in the top two rows especially in cardiovasculature, chromosomal anomalies and body wall defects (Fig. 4).
Supplementary  Fig. 16 shows the geospatial relationships which were derived from spdep::poly2nb and then edited to include all geospatial links. Table 1 gives the results of geospatial regression rstly for a model with additive terms in drug exposure and then a fully interactive model in drug exposure. These are spatial error models and are not spatially lagged. In the additive model series cannabis is independently linked with all eight anomalies particularly cardiovascular (ASD, PDA and tetralogy of Fallot, ToF) and chromosomal (ACD and Downs syndrome), gastroschisis and small intestinal atresia.
In the interactive series of models cannabis is more strongly linked with these same defects. VSD is now positively associated as is diaphragmatic hernia which have both been previously noted to be cannabisassociated [16,17].
A similar exercise is executed for spatially lagged additive (Supplementary Table 11) and interactive (Supplementary Table 12) models with very similar results. In each case spatial error models were superior to combined spatial error and spatial lag (SARAR) models, as judged by the log maximum likelihood values and spatial Hausman tests.
One notes also that in a number of models spatial factors are noted to be highly signi cant. This therefore justi es the use of spatial models and also suggests that spatial factors are signi cant in considering clinical teratological patterns.
Having demonstrated a strong associational relationship between drug exposure and several congenital anomalies the next issue of importance relates to the issue of whether the relationship was causal or not.
Inverse probability weights were generated and used to derive a dataset pseudorandomized for cannabis exposure. Data was processed by robust interactive generalized linear modelling functions. As shown in Table 2 cannabis was signi cantly related to 18 anomalies either alone or in interaction with tobacco and alcohol. It is conceivable that the described relationships were related to some factor other than the measured covariates. E-Values quantitate the degree of association required of some unmeasured confounder with both cannabis exposure and the dependent variables to explain away the described effect. As shown in Table 3 the E-Values were mostly larger than 1.3 and ranged up to 3.8 × 10 30 for geospatial models and up to in nity for mixed effects models, making uncontrolled confounding unlikely.

Discussion
This investigation presents many intriguing ndings. Despite the several technical shortcomings of this dataset it is fascinating for the details and tantalizing clues which have been revealed. Importantly most of its major ndings have been con rmed previously in other locations particularly in Colorado, Hawaii, Canada, and USA and by professional bodies such as AHA, AAP and CDC lending support to the strength of its principal results [2,3,11,13,16,17].
NNSW has higher prevalence rates of the cannabis related anomalies: neural tube defects; small intestinal atresia; body wall defects: exomphalos, gastroschisis, diaphragmatic hernia; the cardiovascular disorders: ASD, VSD, PDA, tetralogy of Fallot, and transposition of the great vessels (TxGrVess); and the genetic disorders: all chromosomal disorders, Downs syndrome, Turners syndrome and trisomy 18. Amongst the defect classes cardiovascular, respiratory, and chromosomal anomalies were elevated.
Some of these associations have been previously reported [3,4,18] and were seen in our unpublished analyses of US data.
QLD Health data showed that the NNSW CI's for CRD's were mostly non-overlapping or were at the extreme end of the QLD CI's. CRD's had higher rate ratios than CNRD's.
Rising rates of cardiovascular, gastrointestinal and respiratory defects, and their rst principal component were associated with falling rates of tobacco and alcohol use but rising cannabis use, just as was found in Colorado and USA [3].
At geospatial and linear regression the cardiovascular defects ASD, VSD, PDA, ToF, TxGrVess; the chromosomal defects ACD, Downs, Turners, Trisomy 13; the body wall defects gastroschisis, exomphalos, diaphragmatic hernia; the GI disorders small intestinal atresia and anal stenosis were all linked with cannabis exposure and for most cannabis exposure was an independent risk factor.
Rising rates of cannabis exposure were more strongly associated with cardiovascular, chromosomal, gastrointestinal and body wall defects than were rising rates of tobacco or alcohol exposure.
Analysis of this dataset by the formal techniques of causal inference analysis including inverse probability weighting and E-Values demonstrated that the described relationships ful l the criteria for causal relationships.
These results show a striking concordance with epidemiological series from elsewhere. ASD, VSD, ToF, obstructive urinary disorders, hydrocephalus, anal anomalies and Downs syndrome were linked with PCE in a large Hawaiian series [4]. VSD has previously been linked with PCE [16]. Neural tube defects were noted to be elevated in a cannabis-related manner in Canada and Hawaii [4,11]. ASD, PDA ACD and Downs were seen to rise in close temporal association with increased cannabis use in Colorado [3].
Exomphalos was implicated in animals [19,20] and in some clinical series including in Queensland [21]. TxGrVess has previously been linked with paternal PCE [22]. Indeed in Canada total CA's were linked with increased cannabis use after controlling for income and sociodemographic variables [2].
Increasing reports from diverse sources indicate that the evidence is building that cannabis has signi cant teratological activities in humans in agreement with animal studies where many severe defects including oedema, exomphalos, phocomelia, spina bi da, myelocoele, exencephaly and foetal loss were documented [19,20]. Concordant reports from Hawaii, Colorado and Canada suggest that the ndings reported herein are indeed valid and are generalizable elsewhere. Given that likely half the NNSW congenital anomalies are reported internally within NSW [7] this suggests that the teratological situation in NNSW is indeed serious. Moreover some of the CA described here, especially chromosomal defects, are heavily therapeutically aborted antenatally again suggesting that the situation may well be much worse than our description suggests. Our analysis strongly implicates cannabis use as a likely underlying factor.
When one also considers the known epigenetic actions of cannabis [2,12,[30][31][32][33] and its associations with developmental neurological dysfunction and autism [34][35][36][37][38]  Moreover since the debate relating to cannabis is typically highly individualistic it seems prudent that medical professional organizations should partner with public health agencies and community groups to enlarge the focus of popular debate from the simply self-referential to a broader multigenerational perspective.
One inevitable conclusion from studies such as this is that access to cannabis should be more highly restricted than at present. Indeed such work calls into question the whole issue of the long term advisability of cannabis medicalization / legalization and the sustainability of such paradigms from a teratological perspective.
The present work has not considered neurological sequalae in the newborn and childhood as has previously been reported to overlap the autistic spectrum disorder and ADHD and thereby potentially play a major role in the modern widespread epidemic of these disorders [34][35][36]39]. When such data is factored into consideration the imperatives for reconsideration and re-evaluation of cannabis legalization overall are largely increased.
Our study has several strengths and limitations. Its strengths include access to whole population data for Queensland and a signi cant portion of the NNSW data. The CA rates and con dence intervals were already provided by QLD Health. The NDSHS is a nationally representative survey conducted every three years and the authoritative source for most Australian drug use data. Our analytical strategy combined CA with drug exposure data which is unusual and useful. We have employed a variety of powerful statistical techniques in this investigation including geospatial analysis, inverse probability weighting, mixed models and E-Values. Study limitations relate mainly to the remote location of the NNSW area close to the Queensland border and the small numbers of some anomalies reported. Losses due to treatment within NSW and to stillbirths and prenatal therapeutic abortion occurring preferentially in CA babies implies that the present ndings are conservative estimates. The very high CA rate reported in Queensland has not been explained despite formal enquiry. The origin of the NNSW denominator gure is unclear. NSW Mothers and Babies reports [7] indicate that during 2008-2015 22,084 babies were born in Northern NSW and 30,848 in the central coast region, totalling 52,932 births. These regions are shown together in our maps. Hence over 11 times the data is available as was used in this analysis if it can be properly collated between the two jurisdictions of NSW and Queensland. This would then facilitate geotemporospatial statistical modelling. This proper collation and assembly of data is a top research priority for future studies. The remote location of NNSW together with its somewhat trans-jurisdictional status has apparently made such a collation di cult in the past.

Conclusions
In conclusion study data indicate that prenatal cannabis exposure is a signi cant and robust covariate of many congenital anomalies in NNSW particularly affecting the cardiovascular, chromosomal, body wall and gastrointestinal systems and is highly signi cant for 10 cannabis-related defects on geospatial analysis. Close concordance between these results and previous reports from Hawaii, Colorado, and Canada and with unpublished USA studies suggest our ndings are reliable and generalizable. Ful llment of the criteria for causal relationships has been demonstrated. Further geospatial epidemiological and basic science research is a priority given cannabis commercialization. Even beyond the obvious jurisdictional health cost-shifting implications careful and thorough further investigation of the teratological pro le of NNSW by coordinated investigations between NSW and Queensland over time to current would appear to be a major international research priority with implications far beyond our shores.

Declarations
Ethics Approval and Consent to Participate Choropleth maps of congenital anomaly class rates across QLD and NNSW. High rates are shown in yellow and low rates in dark blue.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download. SupplementaryMaterial3.docx