|Home | About | Journals | Submit | Contact Us | Français|
To compare gene expression profiles between ankylosing spondylitis (AS) and undifferentiated spondyloarthritis (USpA) patients with inflammatory low back pain.
Peripheral blood mononuclear cells (PBMC) from AS, USpA and healthy subjects were screened with genome-wide microarrays, followed by validation with Real-time PCR.
Microarray profiling and Real-time PCR assays showed that any differences between AS and healthy subjects were only minor. In contrast, 20 genes were strikingly more highly expressed in USpA. The “regulator of G-protein signaling 1” (RGS1) was identified as the most useful biomarker distinguishing especially USpA patients and to a lesser extent AS patients, from control subjects (p= 2.3×10-7 and 6.7×10-3 respectively). All these findings were verified with an independent cohort of patients, which also included rheumatoid arthritis and patients with mechanical low back pain. The Receiver-operator-characteristics (ROC) area under the curve (AUC) values in the first and second cohorts of USpA patients were 0.98 and 0.93 respectively (p=1×10-4). To evaluate the possible derivation of RGS1, we cultured a monocyte-derived cell line with a panel of cytokines and chemokines. RGS1 was significantly induced either by TNF-α or by IL-17.
(1) The PBMC of USpA carry strikingly more highly expressed genes compared to AS or healthy subjects. (2) TNF-α and IL-17 – inducible RGS1 is a potential biomarker for USpA patients and to a lesser extent for AS patients, with inflammatory low back pain.
Spondyloarthritis (SpA) is a family of diseases consisting of the following members: ankylosing spondylitis (AS), undifferentiated spondyloarthritis (USpA), reactive arthritis (ReA), arthritis associated with inflammatory bowel diseases (IBD), and arthritis associated with psoriasis (PsA) (1). The most common are AS and USpA (2, 3). The diagnosis of ReA, PsA and arthritis associated with IBD is frequently made straightforward by the corresponding extra-articular features. However, the diagnosis of AS or USpA can be difficult, especially when the predominant clinical feature is low back pain, and not peripheral extremity inflammation. This is because low back pain alone is a very common primary complaint, contributing to 3% of annual medical visits in the U.S (4). However, only 5% of the chronic back pain seen in general practice is associated with SpA (5). Whereas the back pain associated with AS and USpA is designated as “inflammatory” back pain, non-SpA-associated back pain is commonly designated as “mechanical” back pain (6, 7). Currently, diagnosis of AS and USpA relies on clinical and imaging parameters that, in combination provide an estimate of the degree of probability of whether a patient has SpA (8). The degree of confidence is often not 100%, especially with USpA. This is because the most specific clinical parameter for diagnosing AS in general practice is plain-radiological sacroiliitis. Whereas plain-radiographic sacroiliitis is always positive in AS, it is negative in USpA. Frequently, years of followup are required before some of the USpA patients develop positive sacroiliitis. Both longitudinal and cross-sectional studies indicate that in some patients, USpA is an early form of AS (3, 9). Currently, there is no single blood-derived biomarker which by itself is highly sensitive and specific for distinguishing AS and USpA from patients with “mechanical” back pain (10). Although HLA-B27, TNF-α and perhaps IL-17 are pathogenic factors, as biomarkers when used alone they have either poor sensitivity or poor specificity (11). However, since some of the plain-radiographic sacroiliitis-negative patients show bone edema of the sacroiliac joints by MRI (12), it is possible that such local processes can lead to release of blood-derived protein- or nucleotide- based biomarkers useful in diagnosis of SpA.
In the first part of the study, we used whole-genome 20,589 probe gene expression microarrays each containing 16,283 curated RefSeq genes to screen for genes differentially expressed in AS and USpA relative to healthy subjects. The ultimate purpose of this entire study is to identify a small number of genes whose expression profile might serve as a cost-effective set of surrogate biomarkers for AS and USpA. This particular paper will not follow an alternate strategy which is to identify a diagnostic microarray gene expression pattern.
From the microarray data, we selected 25 most promising individual candidate genes, based on the degree of differential expression and the degree of statistical significance. These were then subjected to a validation test using Real-time PCR, still using the same cohort of AS and USpA patients. At least 6 biomarkers were validated to show potential diagnostic value when assessed by the area under the curve (AUC) of the Receiver-operator-characteristics (ROC). Biomarkers with ROC AUC 0.8 to 1.0 are usually considered as being useful in clinical practice (13). In the third part of our study, we submitted these 6 candidates to an even more stringent level of screening by using an independent cohort of AS and USpA patients as well as healthy subjects. In addition, we also compared the results to those of rheumatoid arthritis (RA) and also patients with mechanical low back pain. From these exhaustive studies, we discovered that the PBMC of USpA contained many more high expression genes compared to AS and healthy subjects. In addition, we identified “regulator of G-protein signaling 1” (RGS1) as the most promising candidate biomarker for USpA and to a lesser extent also for AS patients. In the human peripheral blood mononuclear cells (PBMC), RGS1 has been reported to be constitutively expressed only in monocytes (14). There are as yet no reports on which arthritis-inducing factors might drive the expression of RGS1 in monocytes. To study this, we stimulated in culture a monocyte-derived cell line with a panel of 25 cytokines and chemokines. We discovered that the most effective for generating RGS1 were TNF-α and IL-17. Thus, RGS1 might be a surrogate of the effects of these 2 arthritis-causing cytokines.
The cohorts of subjects providing samples for microarray as well as PCR studies were recruited from Beijing and Guangzhou in China, while the second cohort providing samples for PCR assays was collected from Taichung in Taiwan. All AS and USpA patients also fulfilled the Calin criteria for inflammatory back pain (6). Patients from the first cohort were also selected for their severity of disease activity. Patients from the second cohort were selected without regard to their disease activity or severity. The common selection criterion was clinical predominance of back pain. From Taichung, we also collected a cohort of rheumatoid arthritis (RA) patients, and a cohort of patients with mechanical low back pain. None of the patients in either cohort had psoriasis, inflammatory bowel diseases or history of precedent infections. The study and consent forms have been approved by the respective ethical committees, and signed by all participants.
Peripheral venous blood samples were collected into CPT tubes (BD Diagnostics). PBMC were separated, and RNA extracted with Trizol®. Labeling of RNA was performed using the Illumina TotalPrep® RNA Amplification Kit (Applied Biosystems Ambion, Austin, TX). cRNA was hybridized to Sentrix Human Ref-8_v2 Beadchips followed by staining according to instructions from the manufacturer (Illumina, San Diego, CA).
Data were normalized together by the quantile normalization method (15). Data from AS and USpA were then separately compared to healthy subjects. These microarray results were submitted to 2 consecutive statistical screening processes to eliminate unlikely candidates. In the first screening, we eliminated genes which were less than two-fold different between patients and healthy subjects. In the second screening, the microarray values were log transformed, and genes with low variances among all subjects were again eliminated from subsequent analysis.
Total RNA was reverse-transcribed to cDNA by using QuantiTec Whole Transcriptome Amplification Kit (QIAGEN, Valencia, CA). RPLP0 (ribosomal protein, large, P0, NM_001002) and GAPDH were used as the housekeeping reference genes for PBMC and cell lines respectively. Relative levels of target gene transcripts were assayed in triplicate using quantitative real-time PCR in the 7900HT Fast Real-Time PCR System (Applied Biosystems, Carlsbad, CA). The reactions were performed in 384-well plates (Applied Biosystems) in a volume of 10μl containing SYBR Green PCR Master Mix (Applied Biosystems), 500nM of each primer, and about 10ng of the cDNA template.
With the exception of three pairs of primers, all primers including RGS1 primer pair #1 were purchased from GeneGlobe of QIAGEN. Sequences of primers which were synthesized separately were RPLP0: 5′-CCACGCTGCTGAACATGCT-3′ and 5′-TCGAACACCTGCTGGATGAC-3′ (16); and RGS1 primer pair #2: 5′-TAGTCTTCACAAGCCAGCCA-3′ and 5′GGAAAAACTTCTTGCCAACC-3′ (qPrimerDepot, NIH); RGS-1 primer pair #3: 5′- ATTGAGTTCTGGCTGGCTTG-3′ and 5′-GATTCTCGAGTGCGGAAGTC-3′ (derived from Primer 3 v.0.4.0). Efficiencies of all primers exceeded 90%. Cycling conditions are in the footnote of Table 2. The relative amount of the target gene in each sample was calculated by the 2-ΔΔCt method.
Data analysis utilized the following software: Microsoft Excel, Prism (Graphpad), SPSS and MedCalc (Mariakerke, Belgium), Multiexperiment Viewer (Dana-Farber Cancer Institute, Boston), BRB Array Tools (Biometric Research Branch, NCI), Prediction Analysis with Microarray (PAM), and the 7900HT Fast Real-Time PCR System Software SDS 2.1 (Applied Biosystems).
Expression values of each gene were expressed as mean ± standard deviation. To compare expression levels of the same gene between two groups of subjects, after verifying that the distributions of values were non-parametric, they were submitted to the Mann-Whitney U test. All p values in microarray and in PCR results were corrected by Bonferroni factor for multiple testing (17). The correction factors were 23,000 for microarray and 35 for Real-time PCR data. P values of <0.05 were considered as being statistically significant. In our study, genes which were statistically significant by Mann-Whitney U test also carried false discovery rates of less than 0.05.
U937 and Jurkat cells were cultured with cytokines at the optimum concentrations recommended by the manufacturers (18).
The demographics and clinical parameters of the first cohort of AS and USpA patients are shown in the first cohort columns Table 1. All AS and USpA patients complained of back pain, with only a very few also complaining of pain in the peripheral extremities. All AS and USpA patients were on NSAIDs and some on sulfasalazine. None were on other DMARDs, corticosteroids or biologics. The AS and USpA patients of the first cohort were recruited from two separate institutes. As control subjects, we recruited from the same institutes 13 male and 7 female subjects of mean age 35 ± 9.1 years.
A second cohort of AS and USpA patients was also collected from a third institute to validate the results of the first cohort. The demographics and clinical parameters of the second cohort of patients are also shown in Table 1. The medications they were taking were similar to the first cohort. From the same institute were recruited 26 healthy controls, 13 male and 13 female subjects, with age 27 ± 9 years. We also recruited a cohort of 12 female rheumatoid arthritis (RA) patients of age 43.3 ± 12.6 years. These patients fulfilled the 1988 classification criteria for RA (19). The mean RA DAS28 score was 4.8 ± 0.8 (20); ESR was 27.5 ± 15.7 mm/lst hr; CRP was 0.5 ± 0.4 mg/l. All RA patients were on NSAIDs and methotrexate. None of the patients in the either cohorts were on corticosteroids or biologics. Lastly, we also carefully selected 3 male and 5 female patients diagnosed as having chronic “mechanical low back pain”. Diagnosis was made jointly by an orthopedic and a rheumatology specialist. The mean age was 52 ± 13 years. The diagnosis was based on the clinical parameters as well as either MRI or nerve conduction studies (21).
Hierarchical clustering was used to obtain an overview of gene expression among control, AS and USpA subjects. In the heat map shown in Figure 1, high expression genes were marked in red. The pattern of expression in the USpA group was globally distinct, with a relatively higher number of highly expressed genes compared to either the control or the AS groups.
When the actual microarray values were analyzed, only 4 genes were significantly more highly expressed in AS compared to control subjects: RGS1, NR4A2, HBEGF and SOCS3. The fold-change in AS versus controls was subtle varying from 2.7 to 3.5 fold. Only two genes were under-expressed in AS: DEFA1 and CAMP. The fold-change of mean microarray values as well as statistical significance of these 6 differentially expressed genes is shown in Table 2. In Table 2, empty cells indicate that values corresponding to those cells did not reach statistical significance.
In contrast, in the USpA group, 38 genes were more than 2 fold over-expressed compared to the control group, all with corrected p values <0.05. The 4 genes over-expressed in AS were also over-expressed in USpA. In USpA, as many as 20 genes were more than 6-fold higher than healthy subjects. The list of these 20 most highly differentially expressed USpA genes and their p values are listed in the “microarray column” of Table 2. Nine of those genes appeared to be related to inflammatory processes: CXCL2, CCL20, CCL3L3, CCL3L1, CCL3, IL-1A, IL-1B, IL-8 and IL-6. Out of the 20,589 probes in the microarray, only two genes were exceptional in being over-expressed to more than 3-fold in both AS and USpA. They were RGS1 and NR4A2.
The commonly accepted “gold standard” of gene expression measurements is Real-time PCR. For validation by Real-time PCR, we selected from the microarray data the 6 genes which were differentially expressed in AS, and the 20 most highly expressed genes in USpA. In addition, IL-1RN, IL-6 and SOCS3 were also selected in spite of their lower statistical values, because of their frequently alleged potential in SpA (Table 2).
Genes differentially expressed in USpA but not selected for PCR were: CARD3, MTD118, CXCR4, IRAK3, CTL4, PLAUR, IL1RL2, MAPK6, COX-2, CXCL16, SOD2, CD278, ICAM1, SOCS1 and RGS2. Their microarray values were all less than 4.5 folds higher than control.
Confirming a high degree of accuracy in our microarray procedures, our results showed a high degree of correlation between the mean values of microarray and Real-time PCR data (r=0.75, p<1×10-4). Most of the genes differentially expressed in microarray were also differentially expressed in Real-time PCR. An exception was NR4A2, which was differentially expressed in AS by microarray but not by Real-time PCR. Also confirming the microarray results, other than CCL20, the mean values of all the highly expressed genes were statistically higher in USpA compared to AS (corrected p values = 1.6×10-2 – 7×10-7). Mean values of all Real-time PCR and their degrees of statistical significance are shown in the “Real-time PCR results column” in Table 2. For both USpA and AS, the gene with the highest degree of statistical significance was RGS1 (p= 2.3×10 -7 and 6.7×10-3 respectively). We also ranked the diagnostic potential of each gene by their ROC AUC. RGS1 was also of the highest diagnostic potential with ROC AUC for USpA and AS being 0.99 (95% CI 0.91-0.99) and 0.84 (95% CI 0.69-0.94) respectively (p<1×10-4) (Figure 2).
To parallel what is commonly practiced in patient care in trying to distinguish inflammatory from mechanical low back pain, we combined the USpA and AS into one group designated as “USpA+AS”. The diagnostic potential of all the genes in diagnosing inflammatory low back pain were then ranked by their ROC AUC values. When these calculations were applied to our combined USpA and AS group, only 6 genes showed diagnostic potential. The values of individual subjects are shown in Figure 2. Of these RGS1 still showed the highest promise as a diagnostic biomarker with a very high ROC AUC value of >0.9 (p<1×10-4).
The PCR primers used to measure RGS1 in Figure 2 were RGS1 primer pair #1 which amplified a segment spanning exons 3 and 4 of the RGS1 gene. To ensure that the PCR values in Figure 2 were not artifacts of this particular pair of primers, the samples of AS, USpA and control subjects were assayed again with two different pairs of primers. One pair, RGS1 primer pair #2 amplified a segment overlapping with that of the RGS1 primer pair #1. The third pair, RGS1 primer pair #3 amplified a segment spanning exons 4 and 5 of the RGS1 gene. The correlation coefficient between values derived from the RGS1 primer pairs #1 and #2 was high at 0.88 (p=1×10-6). That for values derived from primer pairs #1 and #3 was also high at 0.89 (p=1×10-6)
RGS1 is a member of the RGS family. To ensure that the results of RGS1 in USpA were unique to the RGS1 member alone, we measured the expression of the following homologues of RGS1 with our USpA and control samples: RGS -2, -3, -4, -5, and -8. Transcripts of RGS -4, -5 and -8 were undetectable. Expression of RGS -2 and -3 in USpA patients were not different from control subjects (p>0.05). These results confirmed the specificity of RGS1 in the RGS family.
The PCR values of all genes shown in Table 2 were higher with the USpA patients compared to the AS patients (p<1×10-4). The relative expression was not due to increased disease activity in USpA relative to AS patients. As a group, these USpA patients were different from the AS patients in having lower values of BASDAI, BASFI and BAS-G (all p<1×10-4), shorter duration of disease (p<3×10-4), better chest expansion and lateral flexion of the lumbar spine, and less radiological sacroiliitis (all p<1×10-4). A matrix correlation analysis showed that there was no positive correlation of RGS1 or any other genes with ESR, CRP, or disease activity indices or peripheral arthritis in either AS or USpA.
We then collected an independent cohort of AS and USpA patients and healthy control subjects from a different geographical area (Table 1). Patient diagnosis was carried out by different rheumatologists. RNA was extracted in a Real-time PCR laboratory instead of the microarray laboratory. Samples were assayed at least 3-6 months after completion of the study of the first cohort. Also included in this cohort was a group of patients with RA, and a carefully selected group of patients with mechanical low back pain (see section on patient demographics).
Samples from this second cohort were tested for Real-time PCR expression of those genes shown in Figure 2, because those genes showed the highest diagnostic potential in the first cohort. To avoid fortuitous omissions, we also tested several genes of only borderline significance in the first cohort. These were CAMP, DEFA1, CXCL2, IL-1α and IL-6. In this second cohort, no differences were detected between either the USpA or the AS groups compared to the control groups for IL-1α, IL-6, CCL3, CAMP and DEFA1. IL-8 and CXCL2 were of poor diagnostic values with ROC AUC <0.6 (Results not shown). As in the first cohort, NR4A2, ATF3 and HBEGF were 6 – 40 folds higher in the USpA group compared to the AS group (P= 5×10-4 – 1×10 -5) (Figure 3).
As in the first cohort, the most useful candidate biomarkers appeared to be RGS1. The mean values of RGS1 in the USpA and AS patients were 20.9 and 11.1 fold higher than those of healthy subjects (p=4×10-6 and 2×10-4, respectively), and were also higher in comparison to both the RA and the mechanical low back pain groups (p= 1×10-2 to 1×10-4). Levels of RGS1 in the latter two groups were not different from those in healthy subjects. Individual PCR values are shown in Figure 3. Even more helpful as a biomarker, for this second cohort, RGS1 appeared to be equally highly expressed in both USpA and AS patients. There were no statistically significant differences between USpA and AS in their mean PCR values or their ROC AUC values. This was in contrast to the first cohort, in which RGS1 was more useful for USpA compared to AS. We also pooled the results of AS and USpA into a single combined group. When this entire group of combined AS and USpA was compared to the mechanical low back pain group, the following biomarker values were observed for RGS1: ROC AUC 0.86, sensitivity 84.2, specificity 87.5, +LR 6.7, PPV 97 and NPV 53.8. Threshold values for these calculations were based on automatic selection by the Medcalc software. The predictive values were based on the assumption that the prevalence of AS and USpA in patients of chronic low back pain was 5.0% (5). The prevalence of SpA among the Chinese population has been reported to be 0.93% (22). When the values of SpA were compared to the general population based on that particular prevalence, the diagnostic values were similar to those for distinguishing from mechanical low back pain.
The RGS1 values shown in Figure 3 were derived from the PCR RGS1 primer pair #1. We also measured the control and the SpA samples with RGS1 primer pair #3. There is a high degree of correlation between results derived from the two different pairs of primers (r=0.85, p=1×10-6).
U937 cells were separately cultured with a panel of 25 cytokines and chemokines, and the expression of RGS1 was assayed 3 and 24 hours later. RGS1 transcription was strikingly induced by culture with TNF-α as well as IL-17 (Figure 4 upper panel). As a control for specificity, we also assayed for NR4A2 in the same samples. In contrast to RGS1, NR4A2 was induced by multiple other cytokines (Figure 4 lower panel). As negative control, we also assayed for RGS1 in Jurkat cells. Similar to a report published previously, RGS1 was undetectable in Jurkat cells (23).
There is a large degree of uncertainty among many clinicians as how to diagnose early AS. That is probably why the diagnosis of AS is frequently delayed for several years (24). Part of the reason is because unlike diseases such as RA, the blood-derived biomarkers currently used for diagnosing AS in clinical practice have low sensitivities and specificities. Numerous papers have been published in attempts to identify more useful blood-derived biomarkers for AS such as MMP3 (10). The present paper is distinct in several regards. First, we based our search on a preliminary screening for differentially expressed genes with genome-wide microarrays. In our analysis, gene selection was driven by statistical significance avoiding the biases of biology-based candidate gene selection approaches. Second, for the microarray study, we tested AS as well as USpA patients within the same study. These patients were recruited by several different investigators, but adhering to the same classification criteria. Third, we verified our microarray results with Real-time PCR. Fourth, we validated the reliability of the most promising candidate RGS1 by comparing its results to homologous members of the same family of genes. Fifth, for the most promising candidate, we validated the results of the first PCR assays by comparing them to two completely different set of PCR primers. The use of these three pairs of primers also suggested that the particular isoform of RGS1 which was expressed in our cohorts was Sp1 (EMBL-EBI Alternative Splicing Database). Sixth, we retested the most promising candidates in an independent cohort of AS and USpA patients and healthy controls, which were recruited by a different clinical investigator from a different geographical location, and also with RNA purified at a different laboratory. Seventh, we compared the AS and USpA results to a cohort of RA patients, and also to a very carefully selected cohort of mechanical low back pain patients.
With each successive experimental step, false candidates were eliminated. Using these procedures, we arrived at three conclusions. The first conclusion is that overall gene expression is higher in USpA compared to AS patients, suggesting that early axial SpA is associated with a more systemic inflammatory process. This is an advantage to clinicians because biomarkers are more in demand in the early stage of SpA rather than the late stage. The second conclusion is that the most promising candidate for both USpA and AS is RGS1. The third conclusion is that, like the other genes observed in this study, RGS1 is probably more highly expressed in USpA compared to AS.
The discovery of RGS1 as an AS and USpA biomarker is unanticipated. Since the expression of RGS1 in our first cohort was much lower in AS compared to USpA, our discovery would have been unlikely without genome-wide screening of both AS and USpA patients. RGS1 is a member of a family of “Regulators of G-protein signaling” (RGS) proteins. These were initially characterized as inhibitors of signaling cascades initiated by G-protein-coupled receptors (GPCRs). RGS proteins are so designated because they can down regulate the signaling initiated by engagement of ligands with GPCRs (25, 26). However, it is now understood that in addition to G proteins RGS proteins also interact with a large network of lipids, ions and other proteins, suggesting that RGS proteins have a diverse array of regulatory functions (27). Among the various RGS proteins, RGS1 is expressed mainly in hematopoetic cells. There is only limited information on the biology RGS1 (28-30). It is not clear from the limited information which particular RGS1 activity is related to SpA disease processes. Nevertheless, although the exact biology of RGS1 is still uncertain, it has already been discovered by random SNP screening to be candidate genes for both Type 1 diabetes and Celiac disease (31, 32), and an independent prognostic marker of disease survival in melanomas (33).
Little is known about the arthritis-related factors which might enhance RGS1 expression (29). We tested a panel of 25 cytokines and chemokines on a monocyte-derived human cell line, and discovered that the two strongest activators of RGS1 expression were TNF-α and IL-17. The kinetics of expression of TNFα and IL-17 were quite different, perhaps reflecting differences in pathways of activation similar to what we observed with HLA-B27 gene (18). Perhaps not fortuitous, TNF-α and IL-17 are currently being considered as key modulators of SpA processes (34, 35).
In conclusion, our genome-wide profiling has identified RGS1 as a candidate biomarker especially with USpA. For any new candidate biomarkers to be used in clinical practice, the candidates will have to undergo the much more vigorous processes of multicenter multi-ethnicity standardization and testing, with well-characterized healthy subjects and also other disease groups (36). Comparison to MRI changes of the sacroiliac joints and the most recent classification criteria for early axial spondyloarthritis will also be essential (37). Our present results serve to contribute a candidate which can be considered for such selection processes. In addition, they might contribute to future discoveries of other more clinical laboratory – friendly and disease-relevant biomarkers, useful in both diagnosis and disease management.
The authors gratefully acknowledge Dr. Chin-Hsiu Liu (Taiwan) for assisting in the Real-time PCR experiments, and Ms Shu-Yuan Chang for collection of clinical data.
Grant supporters: J Gu: grants (30325019 and 30872328) from the National Natural Sciences Foundation of China, clinical subject (2007) of Ministry of public health of China and grants (2002Z3-E4021 and 2005A30801005) from Foundation of Guangzhou and Guangdong province of China; YL Wei and D Yu: Nora Eccles Treadwell Foundation; JCC Wei and MS Jan: National Science Council, Taiwan, NSC 93-2320-B-040-030; F Huang: grant A30771983 from the National Natural Sciences Foundation of China. M.B. Frank:P20 RR016478-077281, P20 RR016478-087197, P20 RR016478-096180, P20 RR020143-047652, U19 AI062629-040006, P20 RR015577-086790 and RR015577-096179 from NIH.
Author Contribution: Dr. Yu had full access to all of the data in the study and takes responsibility for the integrity of the data and the accuracy of the data analysis.Study design. Yu, Frank, Centola.
Acquisition of data: Gu, Huang, Wei, Wei, Jan, Frank
Analysis and interpretation of data: Yu, Frank, Wei.
Manuscript preparation: Yu, Frank, Centola, Wei, Gu.
Statistical analysis: Yu, Frank, Wei.