A Genome-Wide Association Study of spontaneous preterm birth in a European population

Definition and significance of preterm birth

Preterm birth (PTB), defined as birth prior to 37 completed weeks’ gestation, is a major public health problem that affects more than 10% of all births worldwide^1,2. Globally, an estimated 15 million babies are born premature each year^1,2. Despite substantial public health efforts over the past several decades, the U.S. PTB rate remained at 11.72% in 2011^3,4.

PTB is the leading cause of neonatal mortality in non-anomalous newborns^1,3,5,6. It is also associated with a broad spectrum of lifelong morbidity in surviving preterm infants, including neuro-developmental delay, cerebral palsy, blindness, deafness, and chronic lung disease^7–9.

Definition of spontaneous preterm birth

Some PTBs are iatrogenic and can be attributed to obstetric intervention aimed at reducing maternal and/or fetal risk. The remaining PTBs are known as spontaneous PTB (SPTB) and are the focus of research efforts to identify genetic and environmental risk factors. Although SPTB is a pressing health issue, the incomplete understanding of its biology has inhibited development of effective prevention and treatment strategies.

Genetics of spontaneous preterm birth

The etiology of SPTB is complex and multifactorial^9–12, although genetic factors are important contributors^9–14. In addition, PTB prevalence varies among different population groups^15–19. African-American ancestry is consistently associated with an increased PTB risk, even after adjusting for epidemiologic risk factors, such as income, education, lack of prenatal care, and other socioeconomic factors^10,20,21.

A candidate gene approach has identified polymorphisms in genes that encode the progesterone receptor, tumor necrosis factor alpha, interleukins 4, 6, and 10, and mannose-binding lectin^22–33, that are mildly to modestly associated with SPTB. However, results have been inconsistent^29,34–36.

Our genome-wide approach for spontaneous preterm birth

Here, we employed an unbiased, genome-wide approach to search for possible candidate genes associated with SPTB. Genome-wide association studies (GWAS), or whole genome association studies, are a commonly used genetic approach to study a disease or a trait. GWAS compares thousands or even millions of common genetic markers, mainly single nucleotide polymorphisms (SNPs), across individuals with a disease or trait status. There are 1,688 publications identifying 11,299 SNPs that are significant for diseases or traits in 17 different categories^37,38. We obtained the SPTB phenotype and genotype data deposited in the National Center for Biotechnology Information (NCBI) Genotypes and Phenotypes Database (dbGaP)^39,40, study accession phs000103.v1.p1, to perform a GWAS to explore genetic variants associated with SPTB.

Data application and approval

We applied to and received approval from dbGaP^39,40 for access to a dataset for SPTB phenotype and genotype (study accession phs000103.v1.p1). We followed the Data Use Certification Agreement we signed during the application.

To access the data, one must apply and agree to the dbGAP terms of usage. A detailed instructions and procedures for application can be obtained at https://dbgap.ncbi.nlm.nih.gov/.

Data content

This dataset contains participants collected by the Danish National Birth Cohort (DNBC)⁴¹. DNBC is a prospective cohort that enrolled more than 100,000 pregnant women in the first trimester, before most adverse outcomes occurred, and therefore is free from sampling or collection bias. In this dataset, there are approximately 1000 preterm births and 1000 term births as controls. These study subjects were collected from 1997 to 2003. All are singleton gestations. Each birth has records of mother-child pairs. With the exception of 24 children with one or two grandparents from other Nordic countries, all other children in the dataset had parents and all four grandparents born in Denmark.

The case (preterm) group contains births delivered before 37 gestational weeks. The control (term) group contains births delivered at approximately 40 weeks’ gestation. In both preterm and term groups, children born with any recognized congenital or genetic abnormality were excluded. Pregnancies with maternal conditions known to be associated with PTB, iatrogenic or spontaneous (placenta previa, placental abruption, hydramnios, isoimmunization, placental insufficiency, pre-eclampsia/eclampsia), were also excluded by DNBC.

The blood sample (buffy coat) was collected for each mother-child pair. Their whole genomes were genotyped on the Illumina Human660W-Quad_v1_A platform (Illumina, Inc., San Diego, California, USA), which contains more than 500,000 markers. Genotyping was performed by the Johns Hopkins University Center for Inherited Disease Research (Baltimore, Maryland, USA). Further data cleaning and harmonization were done at the GENEVA Coordinating Center at the University of Washington (Seattle, Washington, USA).

Quality control

In this study, we focused only on fetal genomes for further analysis. Individuals with missing genotypes greater than 3% were filtered out. In addition, individuals with a heterozygosity rate deviating more than 3 standard deviations were also excluded from further analysis⁴². After per-individual quality control, we also performed per-SNP filtering. The SNPs that had missing genotypes greater than 3% were excluded⁴². Those having significantly different genotype missing rate between the case and control groups were also eliminated. A conservative cut-off with p < 1×10^-5 was applied⁴². We also excluded SNPs that significantly deviated from Hardy-Weinberg equilibrium in the control group – those with p < 1×10^-5 were filtered out^42,43.

Data analysis

After data quality control, we performed GWAS for the fetal genomes on 22 autosomal chromosomes using the PLINK software package v1.07⁴⁴ (http://pngu.mgh.harvard.edu/purcell/plink/). The level of genome-wide significance was set at 9.18×10^-8, corresponding to Bonferroni correction for 544,675 multiple independent tests.

Quantile-quantile (QQ) plot was made using STATA Statistical Software, Release 12⁴⁵. Manhattan plots were generated using Haploview⁴⁶.

Data overview

We started with a dataset comprised of 1,900 children. There were 31 children having a missing genotype rate greater than 3%. The average heterozygosity was 0.3238, with standard deviation of 0.0059. There were also 31 children having heterozygosity that deviated more than 3 standard deviations. In fact, there was considerable overlap when applying these two filtering criteria (Figure 1) - 26 individuals were identified by both exclusion criteria. During per-individual quality control, a total of 36 individuals were excluded, resulting in 1864 individuals. However, 66 of them had a missing phenotype, yielding a final total of 849 cases and 949 controls (Table 1).

Figure 1. Per-individual quality control.

The X axis is the missing genotype rate for each individual. The Y axis is the heterozygosity rate for each individual. Each dot represents a person. The vertical dash line is the cut-off for per-individual missing rate: 3%. Individuals with missing rate greater than 3% are excluded. The two horizontal dot lines represent the mean heterozygosity ± 3 standard deviations. People with heterozygosity deviate above 3 standard deviations are also excluded. The two criteria overlap largely at the right lower part of the graph.

Among the 560,768 markers, there were 1,933 SNPs that exceeded the missing rate threshold of 3%. There were 367 SNPs that had a significantly different missing rate between the case group and control group (P value < 1×10^-5). Further, 885 SNPs significantly deviated from Hardy-Weinberg equilibrium in the control group (P value < 1×10^-5). These three criteria also identified some overlapping SNPs; a total of 2,670 SNPs were excluded using all criteria. These quality control steps left 558,098 SNPs remaining in the dataset. Of these SNPs, 544,675 of them are located on 22 autosomes, and were included in the analysis (Table 1).

Allelic test

We carried out GWAS for these 1,798 individuals, of which 849 are SPTB cases and 949 are term controls, over 544,675 SNPs on 22 autosomal chromosomes. An allelic test was first carried out, and no SNPs reached genome-wide significance after Bonferroni correction for multiple testing (Manhattan plot (Figure 2) and QQ plot (Figure 3)).

Figure 2. Manhattan plot for allelic test.

The X axis is the position of each SNP grouped by different chromosomes, and presented with different colors. The Y axis is the P value for each test for each SNP, in –log10 scale. The horizontal blue line indicated the threshold of genome-wide significance after Bonferroni correction. The threshold is at 7.04, which corresponds to –log10 (0.05/544,675) for 544,675 independent tests. No SNP reached genome-wide significance.

Figure 3. QQ plot for allelic test.

The X axis indicated the expected P value, in –log10 scale. The Y axis indicated the observed P value from allelic test, also in –log10 scale. The red diagonal line is the line of Y=X, where observed equals expected. The genome-wide significance threshold for –log10 (P) scale is at 7.04, which corresponds to –log10 (0.05/544,675) for 544,675 independent tests. No SNP reached genome-wide significance.

Other genetic models

We then performed GWAS with different genetic models. We tested three classical Mendelian inheritance models here. The recessive model assumes that carrying two variant alleles is required to present a different phenotype; while in the dominant model, one variant allele is sufficient to present a different phenotype as carrying two variant alleles. The additive model assumes the heterozygotes present an intermediate phenotype between the two homozygotes and thus consider the three genotypes separately. The Manhattan plot for the additive model (Figure 4), dominant model (Figure 5), and recessive model (Figure 6) are shown. The dominant or recessive models refer to the action of the minor allele. Within these genetic models, no SNP reached genome-wide significance after Bonferroni correction for multiple testing.

Figure 4. Manhattan plot for additive model.

The X axis is the position of each SNP grouped by different chromosomes, and presented with different colors. The Y axis is the P value for each test for each SNP, in –log10 scale. The horizontal blue line indicated the threshold of genome-wide significance after Bonferroni correction. The threshold is at 7.04, which corresponds to –log10 (0.05/544,675) for 544,675 independent tests. No SNP reached genome-wide significance.

Figure 5. Manhattan plot for dominant model.

This is the GWAS result assuming the dominant action of minor alleles. The Y axis is the P value for each test for each SNP, in –log10 scale. The horizontal blue line indicated the threshold of genome-wide significance after Bonferroni correction. The threshold is at 7.04, which corresponds to –log10 (0.05/544,675) for 544,675 independent tests. No SNP reached genome-wide significance.

Figure 6. Manhattan plot for recessive model.

This is the GWAS result assuming the recessive action of minor alleles. The Y axis is the P value for each test for each SNP, in –log10 scale. The horizontal blue line indicated the threshold of genome-wide significance after Bonferroni correction. The threshold is at 7.04, which corresponds to –log10 (0.05/544,675) for 544,675 independent tests. No SNP reached genome-wide significance.