dataesr/scientific-paragraphs-categorization huggingface dataset:
0017 was only present for one locus (at which significance was achieved) so we examined whether a composite effect of the known cholelithiasis genetic risk variants was detectable in ICP patients by calculating polygenic risk scores (PRS). Of the 29 cholelithiasis loci, genotype data were available in the NIHR-RD dataset for 27. PRS calculated using the OR from the previously published cholelithiasis GWAS demonstrated ICP patients have a significantly higher burden of cholelithiasis genetic risk variants compared with controls (Wilcoxon P = 1. 4 × 10 -05 ). However, when this analysis was repeated with the above loci that showed association in the ICP meta-analysis removed (leaving 18 variants for which data were available), there was no difference in PRS (P = 0. 198) between ICP and controls, despite 92% power to detect a difference at P < 0.
Methods Study cohorts and ethics. An extensive description of the cohorts included in this study and details of the ethics approvals can be found in the Supplementary Note in Methods (Supplementary Information). NIHR-RD and 100KGP whole-genome sequencing variant calling, quality control, ancestry, and relatedness estimation. Detailed methods describing the processing of whole-genome sequence data for the NIHR-RD and 100KGP cohorts Fig. 5  GWAS meta-analysis of ICP reveals association signals affecting the bile acid homeostasis genes CYP7A1 and SULT2A1. a, d LocusZoom  61  plots showing the association signals at CYP7A1 and SULT2A1, respectively. The LocusZoom plots are as per Fig. 3  and association testing was performed as described in Fig. 1. All epigenomic and transcriptomic datasets shown were retrieved from ENCODE  18.
Adult liver ATAC-seq peaks track (green regions) corresponds to regions accessible in at least two out of four ENCODE adult liver ATAC-seq (see "Methods"). 95% PP ICP variants, variants within the 95% genetic credible set in the ICP meta-analysis. LD ICP variants, variants with EUR r 2 > 0. 8 with the meta-analysis lead variant. PP func , functional posterior probability calculated using PAINTOR. CADD, combined annotation-dependent depletion. a Functional fine-mapping with PAINTOR 17 identified a single variant (rs10504255) with a posterior probability over 0. 99 overlapping an active liver enhancer located between the CYP7A1 gene, which encodes the rate-limiting enzyme for bile acid synthesis, and UBXN2B. b, e GTEx 34 [GTEx Analysis Release V8 (dbGaP Accession phs000424. v8. p2)] human liver eQTLs identify the top prioritized variant in the CYP7A1 locus (rs10504255) as an eQTL associated with the eGene UBXN2B (b) and show that the SULT2A1 promoter ICP variant (rs296361) is an eQTL for both SULT2A1 and its downstream lncRNA LINC01595 (e).
NIHR-RD and 100KGP genome-wide association analyses. For genome-wide association analysis (GWAS) in both the NIHR-RD and 100KGP cohorts separately, variants passing the following criteria were retained: MAF ≥ 0. 01, minimum minor allele count (MAC) ≥ 20, missingness <1%, Hardy-Weinberg equilibrium P value >1 × 10 -6 and differential (case/control) missingness P value >1 × 10 -5. The final NIHR-RD and 100KGP datasets comprised 8,337,027 and 9,545,879 variants, respectively. GWAS was performed in each dataset separately using the R package SAIGE  13   weizhouUMICH/SAIGE). SAIGE implements a generalized logistic mixed model to account for population stratification and a saddle-point approximation to control type 1 error rates when case-control ratios are unbalanced. SAIGE was run using the default parameters.
Null logistic mixed models were fitted using the subsets of linkage disequilibrium (LD)-free high-quality common SNVs. Sex and the top ten principal components were used as fixed effects. The genomic inflation factors (lambda), calculated based on the 50th percentile, were 1. 009 and 1. 001 in NIHR-RD and 100KGP, respectively, indicating no evidence of confounding by population structure. The summary statistics from the NIHR-RD were lifted over from build GRCh37 to build GRCh38 of the human genome using CrossMap 59 (version 0. 5. 2) (http://crossmap. sourceforge. net/). In all, 12,756 variants failed to lift over. Due to data sharing restrictions, it was not possible to combine the raw genome-wide genotype data from the NIHR-RD and 100KGP datasets. NIHR-RD and 100KGP meta-analysis. Meta-analysis of the NIHR-RD and 100KGP GWAS summary statistics was performed using METAL (version 2011-03-25) (https://genome.
sph. umich. edu/wiki/METAL)  14. Multi-allelic variants were excluded such that there were 8,291,828 and 9,545,879 variants in NIHR-RD and 100KGP, respectively. Of these, 8,207,819 variants were shared and 8,201,630 variants had matching alleles. Meta-analysis was performed weighting the effect size estimates using the inverse of the standard errors. Variants showing heterogeneity of effect between the two datasets (P < 1 × 10 -5 ) and those in which the minimum/maximum allele frequencies differed by >0. 05 were excluded leaving 8,199,999 variants. A Manhattan plot was produced using the R package qqman (version 0. 1. 8) (https://github. com/stephenturner/qqman)  60 (Supplementary Figure  3a ). Four loci achieved genome-wide significance (P < 5 × 10 -8 ). A QQ plot was generated using the observed and expected -log10 (P values) (Supplementary Fig.
3c ). The genomic inflation factor (lambda), calculated based on the 50th percentile, was 1. 01 indicating no significant population stratification. FinnGen genome-wide association analysis. FinnGen is a public-private partnership combining digital health record data from Finnish health registries with genotyping data from Finnish Biobanks. Release 4 (https://finngen. gitbook. io/ documentation/v/r4/) includes association data at 16,962,023 variants for 2444 endpoints in 176,899 Finnish individuals. A detailed description of the study design and analytical methods are available in the online documentation (URL above). In brief, individuals were genotyped with Illumina and Affymetrix chip arrays. QC was performed to remove samples and variants of poor quality. Genome-wide imputation was performed using reference Finnish whole-genome sequence data.
Disease endpoints were defined using nationwide registries based on ICD (and other) codes. A subset of unrelated individuals of genetically confirmed Finnish ancestry was identified. GWAS was performed using SAIGE (version 0. 35. 8. 8) for all variants with MAC > 5. Sex, age, 10 PCs, and genotyping batch were included as covariates in the analysis. For the ICP GWAS (https://r4. finngen. fi/pheno/O15_ICP), cases (n = 740) were identified as those with the ICD-10 code of O26. 6 ("Liver disorders in pregnancy, childbirth and the puerperium" including "cholestasis (intrahepatic) in pregnancy" and "obstetric cholestasis") or an ICD-9 code of 646. 7 ("Liver and biliary tract disorders in pregnancy"). All females without these ICD codes were included as controls (n = 99,621). Seven loci achieved genome-wide significance (P < 5 × 10 -8 ) (Supplementary Fig.
3b ). The genomic inflation factor (lambda), calculated based on the 50th percentile, was 1. 046 (Supplementary Fig. 3d ). Combined UK (NIHR-RD and 100KGP) and FinnGen meta-analysis. A metaanalysis was subsequently performed of the combined NIHR-RD and 100KGP data ("UK meta-analysis") and the FinnGen summary statistics. This included 154,780 individuals of whom 1138 were cases and 153,642 were controls. Meta-analysis was performed using the same methodology as described above. A total of 7,721,597 variants were shared between the final UK meta-analysis dataset and the FinnGen dataset and 7,715,997 variants had matching alleles. After performing the metaanalysis, variants showing heterogeneity of effect between the two datasets (P < 1 × 10 -5 ) and those in which the minimum/maximum allele frequencies differed by >0.
25 were excluded leaving 7,715,762 variants. Eleven loci achieved genome-wide significance (P < 5 × 10 -8 ) (Fig. 1 ). The genomic inflation factor (lambda), calculated based on the 50th percentile, was 1. 027 (Supplementary Fig. 3e ). Summary statistics for each of the lead variants in each of the separate GWAS and meta-analyses are provided in Supplementary Data 1. LocusZoom. Regional high-resolution association plots showing the LD between markers in the prioritized loci (as per European reference data) were generated using LocusZoom (version 0. 13. 3) (https://my. locuszoom. org/)  61. Credible sets were calculated by LocusZoom using Bayes factors based on the P values (https:// statgen. github. io/gwas-credible-sets/method/locuszoom-credible-sets. pdf). Epistasis. Testing for epistatic interactions was performed across all possible 55 pairs of the 11 prioritized lead risk variants in NIHR-RD and 100KGP separately with PLINK  16  (version 1.
9) (https://www. cog-genomics. org/plink/). The resulting P values were combined using Fisher's method. No significant interactions were identified after correcting for multiple testing based on 55 comparisons (P < 0. 0009). Conditional analysis. Each of the prioritized loci was examined using "conditional and joint analyses" (GCTA-COJO)  15  to assess whether there was more than one independent signal. COJO was run from GCTA (version 1. 93. 2) (https:// cnsgenomics. com/software/gcta/#COJO) using default parameters (examining a 10-Mb window surrounding genome-wide significant variants). The NIHR-RD dataset was utilized to provide the reference LD structure of the variants. Power calculations. Post hoc power calculations were performed to assess the range of odds ratio (OR) and allele frequencies for which associations could be detected at genome-wide significance (P < 5 × 10 -8 ) given the number of cases and controls in the final meta-analysis.
The R package genpwr  62  (version 1. 0. 4) (https:// cran. r-project. org/web/packages/genpwr/index. html) was employed using a logistic model under genetic additivity. Power to detect specific previously described associations was calculated with genpwr using the OR of the described associations, reference European AF and cohort sizes. For cholelithiasis, the OR were taken from the supplementary data of a published meta-analysis  57  (Supplementary Data 1). For rs708686 the OR from the additive model described in the manuscript was utilized. Polygenic risk score calculation. A polygenic risk score (PRS) for cholelithiasis in the NIHR-RD data was calculated using Mangrove (version 1. 21) (https://cran. rproject. org/web/packages/Mangrove/index. html)  63 with the previously published cholelithiasis OR  57  and reference European AF, as described above.
The PRS was calculated assuming additivity within and between loci. Statistical comparison of the PRS in cases and controls was performed using an unpaired two-sample twosided Wilcoxon rank-sum test. Power calculations for the PRS analyses were performed through simulating genotype data for case and control datasets using custom code, written in R (available on request). Genotypes were simulated for each variant in each individual with the probabilities of each genotype status (homozygous reference, heterozygous, and homozygous alternate) determined by the Hardy-Weinberg equilibrium formulae based on the AF and assuming complete linkage equilibrium between the variants. For controls, the reference European AF were used. For cases, the AF of the risk allele were calculated from the reference AF and OR.
A total of 100 repetitions of the simulation process were performed in which genotypes for the target variants were generated for the required number of cases and controls Fig. 6  Analysis human liver eQTLs uncovers TMEM147 and ENPP7 as effector transcripts of ICP susceptibility. a, d LocusZoom 86 plots showing the association signals at TMEM147 and ENPP7, respectively. The LocusZoom plots are as per Fig. 3  and association testing was performed as described in Fig. 1. All epigenomic and transcriptomic datasets shown were retrieved from ENCODE  18. Adult liver ATAC-seq peaks track (green regions) corresponds to regions accessible in at least two out of four ENCODE adult liver ATAC-seq (see "Methods"). 95% PP ICP variants, variants within the 95% genetic credible set in the ICP meta-analysis.
LD ICP variants, variants with EUR r 2 > 0. 8 with the meta-analysis lead variant. PP func , functional posterior probability calculated using PAINTOR  17. CADD, combined annotation-dependent depletion. Overlap of ICP-associated variants in these two loci with adult liver accessible chromatin sites enabled the prioritization of two ICP regulatory variants in intronic regions of GAPDHS (rs4806173) (a) and ENPP7 (rs9916601) (d). b, e The ICP variants prioritized with functional liver chromatin annotations are eQTLs in liver tissue and pinpoint TMEM147 and ENPP7 as effector transcripts in these loci. Violin plots represent the density distribution of the samples in each genotype (n for each genotype is indicated below in blue). Box plots show normalized gene expression in median (white line), first and third quartiles.
Gene-level adjusted P values calculated with FastQTL 86 are shown. Liver eQTL plots were retrieved from GTEx  34  [GTEx Analysis Release V8 (dbGaP Accession phs000424. v8. p2)]. We note that GAPDHS is a testis-specific gene and its interrogation as an eGene for rs4806173 was therefore not possible (b). c The ICP risk allele of rs4806173 (G), which associates with less transcription of TMEM147, is predicted to create a binding site for glucocorticoid receptors, a family of context-dependent transcriptional regulators that can act as either activators or repressors  87   3  and association testing was performed as described in Fig. 1. All epigenomic and transcriptomic datasets shown were retrieved from ENCODE  18. Adult liver ATAC-seq peaks track (green regions) corresponds to regions accessible in at least two out of four ENCODE adult liver ATAC-seq (see "Methods").
8 in European individuals using 1000 Genomes data were identified using LDlinkR (version 1. 1. 2) (https://cran. r-project. org/web/ packages/LDlinkR/vignettes/LDlinkR. html)  64. This resulted in the identification of 871 unique variants. These variants were subsequently analyzed using the Ensembl Variant Effect Predictor (VEP)  65  (version 104) (https://www. ensembl. org/info/ docs/tools/vep/index. html) to identify the 39 unique genes impacted. These genes were interrogated using the GENE2FUNC tool in FUMA (version 1. 3. 6b) (https:// fuma. ctglab. nl/gene2func)  66. 35 genes were successfully mapped and analyzed by GENE2FUNC (ABCB1, ABCB11, ABCB4, ABCG5, ABCG8, ATP4A, BABAM2, C2orf16, CCDC121, CRX, CYP7A1, DHRS9, EIF2B4, ENPP7, GAPDHS, GCKR, GPN1, HNF4A, LINC01595, MRPL33, NSMAF, PPIAP85, PPP4R4, RBKS, SBSN, SDCBP, SERPINA1, SERPINA10, SERPINA2, SERPINA6, SHROOM3, SNX17, SULT2A1, TMEM147, TMEM147-AS1, TPRX2P, UBXN2B, ZNF512, and ZNF513).
Variants within a 100-kb window centering on the lead variant at each of the 11 genome-wide significant loci in the meta-analysis were included. LD matrices of pairwise correlation coefficients were derived using European 1000 Genomes (Phase 3) imputed data  69  , excluding variants with ambiguous alleles (A/T or G/C) where reference/alternate alleles could not be reliably matched. Each variant was intersected with the following functional annotations to generate a PP of being causal (PP func ): PhastCons elements (phastConsElements100way, updated August 5, 2015), ENCODE  18  ATAC-seq replicated peaks generated from the adult liver (see above), and Roadmap E066 (adult liver) ChIP-seq narrow peaks for H3K27ac and H3K4me1 histone modifications. Variants with PP func > 0. 8 were prioritized for further investigation.
For the loci in which PAINTOR did not identify at least one variant with a high posterior probability of being causal (PP func > 0. 8), we then proceeded to intersect the coordinates of all ICP variants in the locus with human adult liver accessible chromatin sites (see details above) using pybedtools version 0. 8. 0  70  (https://daler. github. io/pybedtools/). This approach was taken because none of the remaining loci contained coding variants. To carry out this analysis, we first defined ICP variants as either (a) all variants within the 95% credible set defined in the ICP GWAS meta-analysis by LocusZoom, or (b) all variants in high LD (r 2 > 0. 8) with the meta-analysis lead variant in European (EUR) individuals. High LD variants were identified using LDlink (https://ldlink. nci. nih. gov/?tab=ldproxy).
org/) and pathogenic mutations from the Human Gene Mutation Database (HGMD)  72  (http://www. hgmd. cf. ac. uk/ac/index. php) were used to create genetic tolerance profiles across human meta-domains at amino acid resolution. MetaDome genetic tolerance profiles were previously derived using 56,319 transcripts, 71,419 protein domains, 12,164,292 genetic variants from gnomAD  73  , and 34,076 pathogenic mutations from ClinVar  71  and queried via the online portal at https://stuart. radboudumc. nl/metadome. Interrogation of the coding variant in the GCKR gene was not possible due to a lack of data availability. Finally, we retrieved the CADD scores of pathogenicity (release version 1. 6, GRCh38/hg38) for all coding variants from https://cadd. gs. washington. edu/download. The pathogenicity analysis results are presented in Supplementary Data 3.
Analysis of the survey of regulatory elements data. Survey of Regulatory Elements was previously performed in HepG2 cells  31. In that study, the authors interrogated 5. 9 million human biallelic SNPs for their ability to confer allele-dependent regulatory activity. We retrieved this dataset (SuRE_SNP_table_LP190708. txt. gz) from https://osf. io/pjxm4/ and queried all ICP variants (see details above) at each of the 11 ICP risk loci. Data are presented as non-risk versus risk allele. P values were obtained by a two-sided Wilcoxon rank-sum test, comparing for each SNP the set of SuRE values of all fragments containing the reference allele versus the set of SuRE values of all fragments containing the alternative allele. Following the original report  31  , P < 0. 00173121 were considered as passing a 5% false-discovery rate threshold.
Data mining was performed using Python (version 3. 7) and figures were plotted using seaborn  74  version 0. 11. 2 (https://seaborn. pydata. org/). Transcription factor binding motif analysis. The prioritized ICP-associated variants (see details above) were investigated for their potential effect on TF binding affinity using Motifbreak R version 2. 2. 0 (https://bioconductor. org/packages/release/bioc/ html/motifbreakR. html). To query motifs, we imported MotifDb (https:// bioconductor. org/packages/release/bioc/html/MotifDb. html) (version 1. 30. 0), which contains >4200 TF motifs, including >2800 motifs from studies in human samples and includes the databases HOCOMOCO  75  , HOMER  76  , and ENCODEmotifs  18. For visual representation, motif logos were downloaded from JASPAR 2020  77  (http://jaspar.
genereg. net/). Analysis of human liver eQTLs. Human liver eQTL data shown in this study were retrieved from the Genotype-Tissue Expression (GTEx) Portal 34 on May 12, 2021, which at the time of accession contained data for 208 liver samples. [GTEx Analysis Release version 8 (dbGaP Accession phs000424. v8. p2)]. At each ICP association signal, we queried the lead variant and, if different, the variant that was functionally prioritized. The analysis was restricted to liver eQTLs and eGenes contained in the same TAD as the lead variant (human liver TAD coordinates were retrieved from http://3dgenome. fsm. northwestern. edu/publications. html). Figures were directly downloaded from the GTEx Portal. Analysis of adult liver 3D chromatin conformation. Gene expression and eQTL analyses were restricted to the topologically associating domains (TADs) containing the ICP-associated variants identified in our meta-analysis.
Human liver TAD coordinates (GRCh38/hg38)  78  were retrieved from http://3dgenome. fsm. northwestern. edu/downloads/hg38. TADs. zip. For specific investigation of long-range chromatin interactions in liver tissue, we interrogated two complementary datasets from adult liver: Hi-C  78  and capture Hi-C  79. These datasets were visualized and analyzed using the 3D-genome interaction viewer and database (http://www. 3div. kr/). Visualization of human adult liver datasets. Human adult female liver (right lobe) epigenomic and transcriptomic datasets we downloaded from the ENCODE portal  18  (https://www. encodeproject. org/) in bigwig format and visualized using the UCSC Genome Browser 80 (GRCh38/hg38). All data corresponds to fold change over control for the shown assay. The visualized files had the following identifiers: ENCFF232QBB (ATAC-seq), ENCFF280QYJ (H3K4me3 ChIP-seq) and ENCFF853WPX (H3K27ac ChIP-seq).
Adult liver 3D chromatin interaction datasets were visualized and analyzed using the 3D-genome interaction viewer and database (http://www. 3div. kr/). Tissue-specificity analysis. In order to investigate whether the ICP risk variants at SHROOM3 could be in the vicinity of or interact with liver-specific genes, we used SPECS scores  81  , which were pre-computed for adult liver for all Ensembl (GRCh38. v85) genes using all GTEx 34 samples. In brief, SPECS is a non-parametric tissue-specificity score that is compatible with unequal sample group sizes. SPECS uses all individual data points available and enables the detection of features that are specifically present or absent in one or more tissue types. The SPECS score has been shown to outperform other tissue-specificity scores  81  , including z-score  82  and JSD  83.
Acknowledgements We thank  NIHR BioResource  volunteers for their participation, and gratefully acknowledge  NIHR BioResource Centres, NHS Trusts  and staff for their contribution. We thank the  National Institute for Health Research ,  NHS Blood and Transplant , and  Health Data Research UK  as part of the  Digital Innovation Hub Programme. P. H. D. and C. W. are supported by the  National Institute for Health Research (NIHR) Biomedical Research Centre at Guy's and St Thomas' Foundation Trust and King's College London. A. P. L. is supported by an  NIHR Academic Clinical Lectureship. C. W. is funded by an  NIHR   Senior Investigator award. The views expressed are those of the authors and not necessarily those of the  NHS ,  NIHR , or  Department of Health and Social Care. I. C. is an  Academy of Medical   Sciences Springboard Fellow  ( SBF005\1050 ).
Data availability The genome-wide summary statistics for the final combined meta-analysis carried out in this study have been deposited in the GWAS Catalog under accession code GCP000309, available at https://www. ebi. ac. uk/gwas/studies/GCST90095084. The lead variants' GWAS and GWAS meta-analysis summary statistics generated in this study are provided in Supplementary Data 1. Genotype and phenotype data from the NIHR-RD participants are available from several sources. 4,835 of the NIHR-RD participants were part of the 100,000 Genomes Project-Rare Diseases Pilot. These data can be accessed by application to Genomics England Limited following the procedure outlined at https:// www. genomicsengland. co. uk/about-gecip/joining-researchcommunity/. The genotype and phenotype data from the remaining 7348 NIHR-RD participants can be accessed by application to the NIHR BioResource Data Access Committee at dac@bioresource.