Factors influencing JUUL e-cigarette nicotine vapour-induced reward, withdrawal, pharmacokinetics and brain connectivity in rats: sex matters

Though vaping likely represents a safer alternative to smoking, it is not without risks, many of which are not well understood, especially for vulnerable populations. Here we evaluate the sex- and age-dependent effects of JUUL nicotine vapour in rats. Following passive nicotine vapour exposures (from 59 mg/ml JUUL nicotine pods), rats were evaluated for reward-like behaviour, locomotion, and precipitated withdrawal. Pharmacokinetics of nicotine and its metabolites in brain and plasma and the long-term impact of nicotine vapour exposure on functional magnetic resonance imaging-based brain connectivity were assessed. Adult female rats acquired conditioned place preference (CPP) at a high dose (600 s of exposure) of nicotine vapour while female adolescents, as well as male adults and adolescents did not. Adult and adolescent male rats displayed nicotine vapour-induced precipitated withdrawal and hyperlocomotion, while both adult and adolescent female rats did not. Adult females showed higher venous and arterial plasma and brain nicotine and nicotine metabolite concentrations compared to adult males and adolescent females. Adolescent females showed higher brain nicotine concentration compared to adolescent males. Both network-based statistics and between-component group connectivity analyses uncovered reduced connectivity in nicotine-exposed rats, with a significant group by sex interaction observed in both analyses. The short- and long-term effects of nicotine vapour are affected by sex and age, with distinct behavioural, pharmacokinetic, and altered network connectivity outcomes dependent on these variables.


SUPPLEMENTAL METHODS Functional
The processing of fMRI images was conducted using the open-source RABIES software [1].For both the anatomical and functional images, extra-space around the brain was automatically cropped and temporal spikes were corrected for at each voxel [2].Dummy scans were automatically detected and removed from each EPI.If dummy scans are detected, the median of these volumes provides a volumetric EPI image as reference, given their higher anatomical contrast.Otherwise, a volumetric EPI image was derived using a trimmed mean across the EPI frames, after an initial motion realignment step.Using this volumetric EPI as a target, the head motion parameters were estimated by realigning each EPI frame to the target using a rigid registration.To conduct common space alignment, structural images, were corrected for inhomogeneities, and then registered together to allow the alignment of different MRI acquisitions.This registration was conducted by generating an unbiased data-driven template through the iterative nonlinear registration of each image to the dataset consensus average, where the average gets updated at each iteration to provide an increasingly representative dataset template [3].The finalized template after the last iteration provides a representative alignment of each MRI session to a template that shares the acquisition properties of the dataset, which makes it a stable registration target for cross-subject alignment.After aligning the MRI sessions, this newly generated unbiased template was then itself registered, using a nonlinear registration, to the SIGMA rat brain template [4].To correct for EPI susceptibility distortions, the volumetric EPI was also subjected to inhomogeneity correction, and then registered using a nonlinear registration to the anatomical scan from the same MRI session [5].
Finally, after calculating the transformations required to correct for head motion and susceptibility distortions, transforms were concatenated into a single resampling operation (avoiding multiple resampling) which is applied at each EPI frame, generating the preprocessed EPI timeseries in native space [6].Preprocessed timeseries in common space were also generated by further concatenating the transforms allowing resampling to the reference atlas, at a voxel resolution of 0.3x0.3x0.3mm.
Confound correction was executed on the EPI timeseries resampled to commonspace.Voxelwise linear detrending was first applied to remove first-order drifts and the average image.Motion sources were then automatically removed using a modified version of the ICA-AROMA classifier where classifier parameters and anatomical masks are instead adapted for rodent images [7].The original ICA-AROMA algorithm, which was designed for human data, was adapted to work with rodent data by modifying the hard-coded human priors for anatomical masking and parameter thresholds for component classification.After an initial independent component analysis (ICA) decomposition of the data, four features are extracted from each ICA component spatial map for classification.The component is classified as motion if the CSF or high frequency content fractions are above a given threshold, or if classified by a pre-trained linear classifier combining the brain edge fraction and motion correlation.The CSF mask is inherited from the rodent reference atlas and the edge mask is automatically generated from the brain mask.The threshold for high frequency content was increased as rodents can express higher BOLD frequencies, particularly under medetomidine.The linear classifier was retrained to select new parameters.The automatic process was evaluated by manually classifying motion and network components, derived from a set of scans from the REST-AWK group anesthetized under a medetomidine-isoflurane mixture, and selected parameters to successfully classify clear motion components while avoiding false classification of brain networks.Next, low pass filtering (0.1Hz) and high pass filtering (0.01Hz) was applied [8].Estimated nuisance time-courses during preprocessing were then used for confound regression.More specifically, using ordinary least square regression, the 6 rigid motion parameters, the mean signal from the WM and CSF masks and the global signal were modelled at each voxel and regressed from the data.Before analysis, a spatial Gaussian smoothing filter was applied at 0.3mm full-width at half maximum (FWHM) [8].
Whole-brain connectivity matrices were generated in commonspace for each subject individually using the SIGMA functional template (59 Regions of Interest) by extracting the seed time-course for every parcel and then measuring the cross-correlation (Pearson's r) between every region pair.The correlation values were then re-organized into a whole-brain matrix representing the 'connectivity strength' between every corresponding region pair.

Diffusion
Images were preprocessed using the fMRI Software Library (FSL, v. 6.0.4) and MRtrix (v.3.0.2) [9].Gibbs Ringing Removal [10], followed by PCA denoising [11], was performed first in MRtrix.TOPUP [12,13], followed by EDDY [14], was used to correct for eddy current induced distortions as well as susceptibility-induced distortions.Tractography was then performed using the MRtrix software package.Response functions for single-fibre WM as well as GM and CSF were estimated from the data themselves using an unsupervised method [15].Fibre orientation distribution images were calculated using multi-tissue spherical deconvolution (msmt_csd) followed by images undergoing multi-tissue informed log-domain intensity normalization [16,17].Whole brain tractograms were generated using second-order Integration over Fiber Orientation Distributions (iFOD2) with 10 million streamlines [18], followed by filtering of tractograms [19].
Diffusion data was registered in the same way as the functional data above with one difference: inverse transformation matrices were used to bring the SIGMA atlas regions-ofinterest into diffusion space for final analysis.This was done due to the unique spatial nature of diffusion imaging, and the potential for confounding effects due to resampling and registration of the diffusion data to commonspace.Finally, using the SIGMA ROIs in diffusion space, single subject connectomes weighted by Streamline Count (SC) were produced.

Statistical Analysis
Network-based statistic (NBS) was used to identify significantly different subnetworks (clusters of nodes and edges) between groups [20].Briefly, NBS first identifies edges that surpass a given threshold (suprathreshold links), followed by identification of connected nodes within this subnetwork and finally permutation testing to assign a p-value (controlled for the FWE) to each subnetwork based on its size.Specifically, a test statistic and corresponding p-value is independently computed for each link based on the strength of the pairwise association the link represents.Two tailed t-tests were run to test for main effects of group (p = .05,t-threshold = 3.1) with age and sex introduced to the GLM as a covariate.Statistically significant networks (p < .05)were extracted for analysis of interaction effects.Within these statistically significant networks, the group by age interaction and the group by sex interaction were tested using an f-test (p=.025, f-threshold = 5).

SUPPLEMENTAL DISCUSSION Locomotion
We found that repeated passive nicotine vapour exposure resulted in hyperlocomotion in both adult and adolescent males, but not females.These results match a recent study that compared intravenous nicotine to nicotine vapour; in contrast to i.v.nicotine that produced hyperlocomotion in both male and female rats, passive nicotine vapour only produced hyperlocomotion in males [21].Interestingly, these findings counter several studies that find females are more sensitive to nicotine hyperlocomotion following chronic exposure via other routes, suggesting a unique behavioural result associated with nicotine vapour inhalation [22][23][24].Females did, however, have higher locomotion at baseline; therefore, the lack of hyperlocomotion seen could be a result of a potential ceiling effect.It is likely that the difference in locomotion is a result of differing nicotine levels in the blood as shown by our pharmacokinetic results.Given locomotion's inverted U-shaped curve in response to nicotine [25], the females may be in the latter half of the curve where nicotine is beginning to have a suppressing effect on locomotion.

Weight
Adult females were more sensitive to nicotine's suppression of weight gain compared to adult males as seen previously [26].Also consistent with previous findings, adolescent males were more sensitive to nicotine's suppression of weight gain compared to adolescent females [27].Thus, the effects of nicotine on weight appear to be reliable across exposure routes.

Influence of Estrous and Gonadal Hormones on Behavioural Findings
It is unlikely that the sex differences in nicotine reward, withdrawal, or locomotion observed in the present findings are the result of estrous, as previous studies have found no effect of estrous phase on these behaviours [28][29][30][31][32]. Additionally, the low variability in plasma and brain levels seen within females in our pharmacokinetics results suggest that hormonal differences did not lead to large changes in plasma nicotine or metabolite levels.One study did find that tamoxifen co-administration facilitated CPP in female rats [33], and another found that ovariectomy eliminated nicotine induced CPP [34]; gonadal hormones are thus likely important factors in nicotine reward but the fluctuations during estrous are not enough to cause significant variations in nicotine reward at the sample sizes used in these experiments.It is thus more probable that sex dependent behavioural differences are the result of organizational hormonal effects, as well as overall differences in average gonadal hormone levels.Differences could also result from altered nicotinic acetylcholine receptor (nAChR) response and expression.There are known sex differences in nAChR upregulation in response to chronic nicotine; males show a greater nAChR upregulation which has been seen in mice [35], rats [36], and humans [37].
Additional Pharmacokinetic Considerations A previous study found that adult female rats given repeated i.v.nicotine resulted in a >10fold nicotine plasma concentration compared to males; this difference was attenuated when males were castrated and females were gonadectomized, indicating that gonadal hormones influence nicotine pharmacokinetics [38].In line with this finding, our study found adult females had greater nicotine and nicotine metabolite concentrations than adult males at every time point in both blood plasma and brain supernatant.This was similar in adolescent female brain supernatant but not plasma.Previous studies on vapour pharmacokinetics have been inconsistent; one study found reduced cotinine plasma in females following passive nicotine vapour exposure compared to males [21], while another found increased cotinine in only female adolescents [39].These differences may be a result of the inconsistent power settings used throughout the literature; lower wattage during vapourization results in smaller particle size, leading to significantly increased respirable fraction of aerosol [40].Because larger particles have a lower probability of entering the lungs, it is likely that the size difference in the oropharyngeal cavity between male and female rats results in decreased deposition of large nicotine vapour particles in female or adolescent rats' lungs compared to those of males or adults respectively.It therefore may be important to use low power devices in future experiments when comparing animals that differ in size.Additional variables in the literature such as PG:VG ratio, vapourizer temperature, ohmage, exposure duration, nicotine type (salt or base), and rat strain also likely contribute to varied findings.There were no sex differences in adolescent nicotine plasma levels as has been observed following tobacco smoke exposure or s.c.nicotine injection previously [41].
Hypothalamic Influence on Behaviour Dysfunction in hypothalamic connectivity could also explain the lack of withdrawal observed in the present study as well as sex differences in the long-term effects of adolescent JUUL vapour exposure on sign-tracking in adulthood seen previously in our lab [42].Orexin signaling in the hypothalamic paraventricular nucleus has been shown to be important for the expression of nicotine withdrawal [43].The decreased hypothalamic connectivity seen in females could therefore be inhibiting their ability to show withdrawal-like symptoms.We have also previously found that males but not females demonstrated enhanced levels of signtracking behaviour when previously exposed to chronic nicotine vapour in adolescence [42].Sign-trackers show increased cue-induced activity in the hypothalamus [44]; thus, if this region shows dysconnectivity following chronic nicotine exposure, females could be resilient to adolescent nicotine-induced sign-tracking as a result.

Supplemental Figure 1 :Supplementary Figure 2 :
Arterial plasma nicotine-N'-oxide and norcotinine concentrations 30, and 60 minutes after 10 minutes of JUUL vapour exposure.a) & b) Comparison of adult males to adolescent males.c) & d) Comparison of adult females to adolescent females.Female adults have higher plasma nicotine-N'-oxide at 60 minutes and greater norcotinine at 30 minutes compared to adolescents.e) & f)Comparison of adult male to adult females.Adult females have higher nicotine-N'-oxide at 60 minutes and greater norcotinine at 30 and 60 minutes compared to adult males.g) & h) Comparison of adolescent male to adolescent female plasma concentrations.*p<0.05adult versus adolescent or male versus female.N=7 per group and timepoint.Data presented as mean ± SEM.Common resting-state networks identified using spatial group-level ICA.The images are shown in an anatomical view overlaid on a commonspace anatomical template, resampled to the EPI's dimensions and the networks are scaled in pseudo-z-scores, i.e., how many standard deviations away from the background noise.

Table 2 :
All connections identified by NBS statistics to have statistically significant interaction Group by Sex effect (p < 0.001, 5 edges, 6 nodes).