This study combined metagenomic sequencing, metabolite measurements, biological-age markers and computational analyses across multiple human cohorts to investigate genetic variation in the gut microbiome and its relationships with host traits. The analyses included Dutch, Tanzanian, inflammatory bowel disease, obesity and functional-genomics cohorts, covering microbial single-nucleotide variants (SNVs), structural variants, population structure, metabolites, telomere length, immune-cell markers and host age. The researchers also modelled protein structures, examined a conserved glutathione-transport gene cluster in Faecalibacterium, identified genes potentially involved in nicotine degradation and developed a microbial SNV-based age-prediction model.
Cohorts
Dutch Microbiome Project
The Dutch Microbiome Project (DMP) includes 8,719 individuals and forms part of the Lifelines study, a multidisciplinary, prospective, population-based cohort that uses a three-generation design to examine health and health-related behaviours in 167,729 people living in the northern Netherlands. Lifelines uses biomedical, sociodemographic, behavioural, physical and psychological assessments, with a particular focus on multimorbidity and complex genetics48. Details of DMP microbiome-data generation are reported by Gacesa et al.12. Briefly, fresh-frozen faecal samples were collected, and microbial DNA was isolated using the QIAamp Fast DNA Stool Mini Kit (FSK, Qiagen) on the QIAcube automated sample-preparation system (Qiagen). Metagenomic sequencing was performed by Novogene, China, using the Illumina HiSeq 2000.
Lifelines-DEEP cohort
Lifelines-DEEP (LLD) is a Lifelines sub-cohort established to investigate risk factors for complex diseases. It includes 1,539 participants with multiple layers of omics data collected at baseline in 201313. High-quality metagenomic data were available for 1,135 participants: 474 men and 661 women. Their mean age was 45.04 years (18–81 years, s.d. = 13.60 years), and their mean BMI was 25.26 kg m−2 (16.67–48.56 kg m−2, s.d. = 4.18 kg m−2). Microbiome-data generation is described by Zhernakova et al.49. Fresh-frozen faecal samples were collected, DNA was isolated with the AllPrep DNA/RNA Mini Kit (APK, Qiagen, 80204), and sequencing was performed on the Illumina HiSeq platform at the Broad Institute in Boston. All 1,135 metagenomic samples passed quality control.
Follow-up APK samples were collected from 338 LLD participants: 150 men and 188 women. The interval between sampling points ranged from 3.33 to 3.92 years (mean = 3.53 years, s.d. = 0.12 years). At baseline, the mean age was 48.2 years (18–80 years, s.d. = 11.7 years) and the mean BMI was 25.4 kg m−2 (17.6–43.3 kg m−2, s.d. = 4.08 kg m−2). At follow-up, the mean age was 51.7 years (22–84 years, s.d. = 11.7 years) and the mean BMI was 25.6 kg m−2 (16.1–37.6 kg m−2, s.d. = 4.0 kg m−2). Two aliquots of each stool sample were processed for DNA isolation and library preparation with the FSK and APK kits. Whole-genome sequencing (WGS) was performed on the Illumina HiSeq 2000.
500-Functional Genomics cohort
The 500-Functional Genomics cohort (500FG) is part of the Human Functional Genomics Project (HFGP) and includes 534 individuals. Among 517 500FG-FSK participants, the mean age was 28.47 years (18–75 years, s.d. = 13.69 years), and the mean BMI was 22.75 kg m−2 (15.09–35.43 kg m−2, s.d. = 2.85 kg m−2). Metagenomic-data generation was described previously14,15. DNA was isolated from 517 and 471 faecal samples using the FSK and APK kits, respectively, and libraries were sequenced on the Illumina HiSeq 2000.
300-Obesity cohort
The 300-Obesity cohort (300OB) was established by Radboud University Medical Center in Nijmegen, the Netherlands50. It includes 298 older individuals with metagenomic data: 132 women and 166 men. Participants had a BMI greater than 27 kg m−2. Their mean age was 67.07 years (54–81 years, s.d. = 5.39 years), and their mean BMI was 30.73 kg m−2 (26.28–45.47 kg m−2, s.d. = 3.48 kg m−2). Participants were enrolled between 2014 and 2016. Faecal metagenomic sequencing used a protocol and analysis pipeline similar to that used for LLD.
1000 Inflammatory Bowel Disease cohort
The 1000 Inflammatory Bowel Disease cohort (1000IBD) consists of patients recruited at the specialist IBD outpatient clinic of the University Medical Center Groningen (UMCG), the Netherlands16,51. IBD diagnoses were based on accepted radiological, endoscopic and histopathological assessments. The cohort was 60.70% female and 39.30% male. The mean age was 42.61 years (18–82 years, s.d. = 14.69 years), and the mean BMI was 25.44 kg m−2 (15–49.4 kg m−2, s.d. = 5.05 kg m−2).
300-Tanzania Functional Genomics cohort
For analyses involving non-European participants, the study used 300TZFG, an HFGP population cohort of 323 individuals from rural and urban areas of Tanzania’s northern Kilimanjaro region. Metagenomic-data generation was described previously18. Bacterial DNA was isolated with the APK kit, and libraries were sequenced on the Illumina NovaSeq 6000. Among participants with metagenomic data, the mean age was 33.08 years (17.76–66.43 years, s.d. = 11.25 years), and the mean BMI was 24.98 kg m−2 (16.41–41.91 kg m−2, s.d. = 4.93 kg m−2).
Metabolite measurements
Bile acids
Fifteen bile acids and C4 concentrations were measured in fasting plasma by liquid chromatography–tandem mass spectrometry (LC–MS/MS), as previously described52,53. For preparation, 25 µl of plasma was mixed with 250 µl of internal-standard solution to precipitate proteins. Samples were centrifuged at 15,800g, and the supernatant was transferred to a clean glass tube. The solvent was evaporated under nitrogen at 40 °C, and samples were reconstituted in 200 µl of 50% methanol in water.
Quantification used a Nexera X2 ultra-high-performance liquid chromatography system (Shimadzu) coupled to a SCIEX QTRAP 4500 MD triple-quadrupole mass spectrometer. The system was controlled with Analyst MD 1.6.2 software. Bile acids were separated on an ACQUITY UPLC BEH C18 column (1.7 µm, 2.1 mm × 100 mm) with an ACQUITY UPLC BEH C18 VanGuard pre-column (1.7 µm, 2.1 mm × 5 mm; Waters). Separation took 28 minutes using 10 mM ammonium acetate in 20% acetonitrile as mobile phase A and 10 mM ammonium acetate in 80% acetonitrile as mobile phase B, at a flow rate of 0.4 ml min−1.
Proportions of the 15 bile acids, identified with the suffix “_p”, were calculated by dividing each concentration by the total bile-acid concentration. Eight additional totals or ratios were calculated:
- Total_BAs: the sum of all bile-acid concentrations.
- Total_primary_BAs: the sum of all primary bile-acid concentrations.
- Total_secondary_BAs: the sum of all secondary bile-acid concentrations.
- Secondary_primary_ratio: Total_secondary_BAs divided by Total_primary_BAs.
- CA_CDCA_ratio: (CA + TCA + GCA) divided by (CDCA + TCDCA + GCDCA).
- The unconjugated-to-conjugated bile-acid ratio: (CA + CDCA + DCA + LCA) divided by (TCA + GCA + TCDCA + GCDCA + TDCA + GDCA + TLCA + GLCA).
- CA_dehydro_deconju_ratio: (DCA + TDCA + GDCA) divided by (CA + TCA + GCA).
- Taurine_glycine_ratio: (TCA + TCDCA + TDCA + TLCA) divided by (GCA + GCDCA + GDCA + GLCA).
Short-chain fatty acids
Short-chain fatty acids in stool were measured by gas chromatography–mass spectrometry, as previously described54,55.
Trimethylamine N-oxide and precursors
Trimethylamine N-oxide and its precursors betaine, carnitine, deoxycarnitine and choline were measured by LC–MS/MS, as previously described56,57. An online sample-cleanup step used an analytical column containing a C18 stationary phase before the sample entered the mass spectrometer. The stationary phase did not retain the analytes but did retain substantial matrix interferents, including (phospho-)lipids58.
Untargeted blood metabolome
Untargeted metabolic profiling of fasting plasma was performed at General Metabolics using flow-injection time-of-flight mass spectrometry on an Agilent 6550Q-TOF system, as previously described59,60,61.
Biological-age measurements
Telomere length
Telomere length was measured with an automated multicolour flow fluorescence in situ hybridization (Flow-FISH) method, as previously described62,63. White blood cells were obtained by lysing erythrocytes in whole blood with NH4Cl. The cells were combined with bovine thymocytes of known telomere length, which served as an internal reference. After denaturation in formamide at 87 °C, the cells were hybridized with a fluorescein-labelled (CCCTAA)3 peptide-nucleic-acid probe specific to telomere sequences and counterstained with LDS751 DNA dye.
Flow cytometry measured fluorescence in granulocytes, total lymphocytes and lymphocyte subsets identified with antibodies against CD20, CD45RA and CD57. Median telomere length was determined from duplicate measurements by comparing fluorescence with internal-control and unstained cells.
sjTREC
Signal-joint T-cell receptor excision circles (sjTRECs), produced during T-cell development, were quantified by quantitative PCR (qPCR), as previously described63. DNA was extracted from whole blood by proteinase K digestion, phenol extraction and isopropanol precipitation. TaqMan qPCR was performed on a ViiA 7 system (Life Technologies), using 75–100 ng of DNA in a 12 µl reaction mix. The reaction contained 18 µM sjTREC primers (5′-TCGTGAGAACGGTGAATGAAG-3′ and 5′-CCATGCTGACACCTCTGGTT-3′), albumin primers (5′-TGAACAGGCGACCATGCTT-3′ and 5′-CTCTCCTTCTCAGAAAGTGTGCATAT-3′), and 5 µM probes for sjTREC (5′-(FAM) CACGGTGATGCATAGGCACCTGC-3′ (TAMRA)) and albumin (5′-(FAM) TGCTGAAACATTCACCTTCCATGCAGA-3′ (TAMRA)).
The PCR programme was 50 °C for 2 minutes, 95 °C for 15 minutes, followed by 45 cycles of 95 °C for 15 seconds and 60 °C for 30 seconds. Reactions were run in duplicate. Probe-cleavage fluorescence was measured at the threshold cycle (CT); valid results had a duplicate s.d. of no more than 1. Normalized sjTREC expression was calculated from the difference between the mean albumin and sjTREC CT values.
Metagenomic processing
Sequencing quality control
Host-contaminated and low-quality reads were removed with KneadData v0.7.4, Bowtie2 v2.3.4.364 and Trimmomatic v0.3965. Raw reads were first aligned to the human reference genome GRCh37–hg19 to remove human DNA. Adapter sequences and low-quality reads were then removed with Trimmomatic’s default settings: SLIDINGWINDOW:4:20 and MINLEN:50.
Reference-genome abundance
Reference-genome abundance in metagenomic samples was estimated with CoverM v0.7.066, using the same reference genomes as those used for metagenomic SNV detection.
Structural variants
Structural variants were defined as highly variable genomic fragments longer than 1 kb occurring within a species across microbial genomes. Gut-microbiome structural variants were detected from cleaned reads with SGVFinder v2.04, again using the reference genomes employed for SNV detection.
SNVs
Metagenomic SNVs were characterized with GT-Pro v1.0.1 and the pre-constructed sck-mer database11. A catalogue of 104 million biallelic SNVs in core genomic regions of 909 gut microbial species was first constructed from Unified Human Gastrointestinal Genome (UHGG) assemblies. Species-specific k-mers were used to tag SNV genotypes and create a k-mer-based metagenomic SNV panel. GT-Pro genotyped gut microbiota by searching for and matching genotype k-mers in metagenomic sequencing data. The output was converted to SNV-dosage format, in which the genotype represented the frequency of the alternative allele.
SNV quality control
For each SNV across n samples, Ri represented the proportion of the reference allele in sample i, and Ai represented the proportion of the alternative allele. These proportions generally sum to 1 when only two alleles are considered. Reference allele frequency (RAF), alternative allele frequency (AAF) and minor allele frequency (MAF) were calculated as:
RAF = ∑i=1NRi / N
AAF = ∑i=1NAi / N, or AAF = 1 − RAF
MAF = min(RAF, AAF)
Quality control was performed separately in each cohort. SNVs with a missing rate above 90% or a MAF below 0.05 were removed, producing cohort-specific common SNVs. Species with fewer than 1,000 common SNVs were also removed. Genotype-data manipulation and quality control used PLINK2 v2.0.0-a.6.9LM67.
Gene and SNV annotation
Genes in gut-microbial reference genomes were predicted and annotated with Bakta v1.10.1 using database v5.168. Functional effects of common SNVs were annotated with SnpEff v5.1f21.
Population structure and association analyses
Genetic relatedness
For each species, between-sample genetic dissimilarity Di,j was calculated as a Manhattan distance scaled by the number of shared SNVs in each sample pair, using the parallelDist() function in the R package parallelDist v0.2.6. Pairs sharing fewer than 1,000 SNVs were excluded and their dissimilarities marked as missing. To retain the largest possible sample without missing values, the sample with the most missing values was removed from the matrix in each round until no missing values remained.
Genetic relatedness was then defined as Ri,j = 1 − Di,j. The matrix MR containing these values represented population genetic structure. Principal-coordinate analysis was performed on the dissimilarity matrices MD with the pcoa() function in the R package ape v5.7-1, using correction = “cailliez”.
Linkage disequilibrium
To summarize linkage disequilibrium (LD) in each species, 1,000 SNVs were randomly extracted from the cleaned genotype dataset with PLINK2 v2.0.0-a.6.9LM using -thin-count 1000. LD was calculated as the squared correlation for every pair of extracted SNVs with -r2 inter-chr yes-really -ld-window-r2 0. Genome-wide LD structures were visualized with LDBlockShow v1.4069, using -SeleVar 2 -NoShowLDist 100000000 -OutPdf -OutPng.
Quasi-phaseability
Multiple strains of the same species within a host can affect population-structure estimates and reduce the power of association analyses. A sample was considered quasi-phaseable when more than 80% of its SNVs had the same allele in more than 90% of reads. This indicates that the sample was probably dominated by one strain of the species. The concept and criterion were described previously38,70.
SNV-based metagenome-wide association study
Phenotypes from the different cohorts and datasets were manually integrated and curated to ensure consistent units. Each phenotype was analysed with a linear mixed model (LMM) containing fixed effects for the corresponding species abundance, age, sex and sequencing depth, together with a random effect representing within-species population structure. The random-effect variance–covariance matrix was the sample-by-sample kinship matrix MR. When age or sex was the analysed trait, that variable was omitted from the covariates.
The LMM analyses used EMMAX beta-07Mar201071 with the parameters emmax -d 10 -t genotype.file -p phenotype.file -c covariate.file -k kinship.file -o output -Z. Analyses were performed separately for each cohort and species, using only quasi-phaseable samples for the relevant species.
For SNV–phenotype pairs available in at least two independent cohorts—DMP, LLD1, 500FG, 300OB, IBD and 300TZFG—meta-analysis was conducted with Metal v2020-05-0572. LLD2 and 500FG-APK were excluded from this meta-analysis. Metal used its default settings, weighting cohort-based P values by sample size. Robustness was assessed with a permutational LMM using 999 permutations and the Monte Carlo procedure within each cohort. Regional association plots showing LD around associated gene regions were generated with LDBlockShow v1.4069.
Significance thresholds
The study-wide significance threshold was based on 922,497,753 tested phenotype–SNV pairs and a Bonferroni correction: 0.05 divided by the number of pairs = 5.42 × 10−11. Suggestive significance was likewise calculated with a Bonferroni correction based on the number of SNVs tested for each species and phenotype. This suggestive threshold did not account for the number of species and phenotypes in the study.
LD-based grouping
Because SNVs can be correlated through LD, related SNVs were grouped to estimate the number of independent association signals. SNVs with between-SNV r2 above 0.01 within a 10,000 kbp window were grouped with PLINK2 v2.0.0-a.6.9LM67. The most significant SNV in each LD block was selected as the representative or tag of an independent signal. PLINK2 used the study-wide threshold of 5.42 × 10−11 and species-specific significance thresholds as its primary and secondary cut-offs. Local LD was visualized with LDBlockShow v1.4069.
Power analysis
A post-hoc approximate power analysis was performed for tag SNVs associated with significant findings using the R package genpwr v1.0.4. The proportion of phenotypic variance explained (R2) was estimated from the observed effect size, standard error, sample size and MAF with pheno_var_explained(). Power was calculated with power.calc.linear() under an additive genetic model and a study-wide α of 5.42 × 10−11. Power was also estimated at fixed sample sizes of n = 100, 500, 1,000, 5,000 and 10,000.
Protein-structure analysis
Because AlphaFold3 supports protein-complex modelling, the three-dimensional structures of the OPP complex were modelled with AlphaFold373 and visualized in ChimeraX v1.10.174 and PyMOL v2.5.7. AlphaFold3 server terms of use prohibit using its output in docking or screening tools (https://alphafoldserver.com/about). The researchers therefore used local AlphaFold2 for monomer prediction, docking and downstream analyses.
To examine interactions between OppA/SarA_2 and glutathione, wild-type and mutant OppA/SarA_2 sequences were modelled with AlphaFold2 v2.3.1 under monomer settings. Default reference databases and a template cut-off date of 14 May 2020 were used. Five models were generated for each experiment and refined with GPU-based AMBER relaxation. The model with the best pLDDT score, with a mean of at least 90, was selected after excluding low-confidence N-terminal residues.
Glutathione was downloaded from PubChem and converted to PDB format (https://pubchem.ncbi.nlm.nih.gov/). P2Rank v2.576 predicted ligand-binding pockets. Residues with ligand-binding probabilities above 0.6 in top-scoring pockets with probabilities above 0.85 were used as active sites in HADDOCK v2.477,78. HADDOCK generated 200 complexes per variant. The four highest-ranked poses from the best-scoring cluster were selected and energy-minimized in water. PRODIGY-LIG v1.1.479,80 estimated binding free energies (ΔG), and PLIP v202581 identified non-covalent interactions.
For benchmarking, a literature-derived82,83 tripeptide panel was docked to the same OppA binding site with HPEPDOCK v2.084. The panel represented basic peptides KAK, KWK and KKK; neutral peptides AAA, GFG and GAG; and acidic peptides EEE and DED. The same active residues used for glutathione docking were applied. Complexes were refined by water minimization in HADDOCK v2.4, and binding free energies were estimated with PRODIGY-LIG v1.1.4.
Structural annotation of SNVs
Protein-structure analysis used an internal pipeline (https://github.com/GArcho44/snvs-annotation-pipeline). Only missense variants in coding regions were retained. Bakta annotations68 and UniProt information (https://www.uniprot.org/)85,86 linked genes to protein sequences in UniRef100 and to structures in the AlphaFold Database (AlphaFoldDB v2022)27.
For genes without direct UniRef100 matches, coding sequences were extracted with GffRead v0.12.787 and translated with EMBOSS transeq v6.6.088. Protein sequences were searched against AlphaFoldDB and the ESM Metagenomic Atlas (1 November 2022; https://esmatlas.com/)26 using DIAMOND BLASTP v2.0.1589. Only full-length matches with at least 98% identity were retained. When multiple matches were available, the match with the highest identity was selected, with AlphaFoldDB preferred in ties. UniProtKB or MGnify accessions were used to retrieve high-confidence models with mean pLDDT values of at least 80 through custom application programming interfaces.
Relative solvent accessibility was calculated with DSSP v3.0.090,91, with surface exposure normalized to empirical per-residue maxima. Ligand-binding pockets were detected with fpocket v4.2.192 and rescored with fpocketPRANK. Binding residues were defined from pockets with druggability scores above 0.5. Each residue’s distance to ligand was calculated as the minimum Euclidean distance to the nearest pocket centre.
Protein–protein interfaces were predicted with PeSTo, model i_v4_193. Residues with interface-confidence scores above 0.5 were clustered by spatial proximity using a Cα–Cα distance below 10 Å. Only high-confidence clusters with an average score above 0.8, pLDDT above 0.7 and at least six residues were retained. Each residue’s distance to the protein–protein interface was calculated as its Euclidean distance to the nearest predicted interface residue.
The three-dimensional structure of a marR homologue dimer was modelled with AlphaFold373 and visualized with ChimeraX v1.10.174. Protein–protein interface analysis, including in silico alanine scanning, used the mCSM-PPI2 web server94.
Glutathione transport gene cluster
The original gene annotation used Prokka v1.13.395. An age-associated signal was identified in a gene annotated as sarA (UniProtKB: P31306), located next to oppB, oppC, oppD and oppF. sarA encodes a periplasmic oligopeptide-binding protein, SarA, and contains an oppA family domain. According to previous studies96,97, OppABCDF is a typical oligopeptide permease, or ABC transport, complex. The five genes could therefore form a cluster whose products assemble into an oligopeptide ABC-transport system.
To test whether this organization was conserved across the Faecalibacterium prausnitzii taxonomic clade, the researchers searched 17,496 publicly available and in-house F. prausnitzii assemblies with cblaster v1.3.2098. The assemblies came from UHGG v2.0.299, SPIRE100 and in-house F. prausnitzii metagenome-assembled and isolate genomes6. All five genes co-occurred in 2,185 genomes, with a maximum distance of less than 20 kb between any two hits in a cluster.
The genomes containing the complete OPP operon were assigned to Faecalibacterium butyricigenerans (223), Faecalibacterium duncaniae (3), Faecalibacterium longum (1,568), F. prausnitzii _I (372), Faecalibacterium sp900539945 (1), Faecalibacterium sp934514285 (15) and three unclassified genomes. The results indicated that sarA, oppB, oppC, oppD and oppF form a conserved gene cluster or operon within the F. prausnitzii genome.
AlphaFold3 was also used to model the protein complex. The predicted structure resembled experimentally determined ABC-transporter structures from Mycobacterium tuberculosis (PDB: 8JA7)24 and Escherichia coli MetNI (PDB: 6CVL)25. The analysis therefore indicated that SarA, through its OppA domain, can function as OppA and assemble with the other four subunits into an OPP complex.
To investigate a possible ligand, gene functions were annotated with Bakta v1.11.068. One of the five genes was annotated as an ABC-type glutathione-transport-system ATPase component (gsiA). The protein sequences of SarA/OppA, OppB, OppC, OppD and OppF were compared by BLASTP v2.16.0+ with 3,239 existing protein sequences from ABC-type glutathione transport systems. All five proteins matched known sequences from this system at E-values below 1.00 × 10−10. The OPP complex may therefore function as a glutathione transport system.
Nicotine-degradation genes
To identify genes potentially involved in nicotine degradation, the researchers obtained sequences for 92 proteins involved in two nicotine-degradation pathways from the Kyoto Encyclopedia of Genes and Genomes (KEGG; entries M00810 and M00811). These sequences were searched against genes in 909 reference genomes with tblastn v2.14.1-gompi-2023a101. In total, 2,434 genes similar to known nicotine-degradation genes were identified at E-values below 1.00 × 10−10.
Microbial SNV-based age prediction
A host-age prediction model was built with a light gradient-boosting machine implemented in the R package lightgbm v3.3.5. Generalizability was assessed with nested cross-validation in tidymodels v1.1.0, using a five-times-fivefold outer split and a fivefold inner split.
To reduce information leakage between individuals with closely related microbiomes, samples from the same lineage group were assigned to the same split—either the internal training set or the internal test set—using group_vfold_cv(). Lineage groups were defined from all SNVs used to train the age model, pooled across species. Pairwise scaled Manhattan distances D were calculated and converted to kinship values K = 1 − D. A neighbour-joining tree was constructed from D and rooted to orient the cut. Lineage groups were assigned using a single tree-height cut of h = (1 − c)/2, with a kinship cut-off of c = 0.9. Thus, h = 0.05, equivalent to a distance cut-off of 0.1. Each descendant clade crossing the cut was treated as one lineage group.
All cross-validation runs used predetermined hyperparameters: trees = 2000, tree_depth = 15, min_n = 30, sample_size = 0.6, mtry = round(0.2 × total_feature_number), and learning_rate = 0.01. Feature importance was estimated with Shapley additive explanations, which assign an importance value to each feature for an individual prediction. After 5 × 5 rounds of fivefold cross-validation, the median importance of every SNV was calculated. SNVs were ordered from highest to lowest median importance, and features were cumulatively included until they accounted for 95% of the total importance. These features were used in the final model.
Because lineage groups were defined from SNVs pooled across species, the procedure captured overall microbiome kinship but might not completely prevent information leakage for individual species whose species-specific kinship differs from the overall pattern.
Partial correlations
Partial correlations were calculated with the pcor() function in the R package ppcor v1.1. Associations between gut-microbial SNV age and other biological-age measures, and between gut-microbial SNV age and intrinsic host phenotypes, used partial Spearman correlations adjusted for chronological age. Associations between gut-microbial SNV age and the untargeted blood metabolome used partial Pearson correlations adjusted for chronological age, sex and blood creatinine concentration.
Ethical approval
The Lifelines study was approved by the UMCG Ethics Committee (METc 2007/152). All participants provided informed consent before enrolment, and DMP participants or their legal representatives provided additional written consent for children younger than 18 years.
The LLD study was approved by the UMCG Institutional Ethics Review Board (ref. M12.113965), the 300OB study by the Institutional Review Board of CMO Regio Arnhem-Nijmegen (nr. 46846.091.13), and the 500FG study by the Ethical Committee of Radboud University Nijmegen (NL42561.091.12, 2012/550). Volunteer inclusion and experiments followed the principles of the Declaration of Helsinki, and all volunteers provided written informed consent before material was collected.
The 1000IBD study was approved by the UMCG Institutional Review Board in Groningen (ref. M12.1139652008.338). The 300TZFG study was approved by the Ethics Committee of KCMC University (KURERC; No. 936) and Tanzania’s National Institute for Medical Research (NIMR/HQ/R.8a/Vol. IX/2290).
Reporting summary
Additional information about the research design is available in the Nature Portfolio Reporting Summary linked to the article.
Source: www.nature.com