Chromosomal Region 11p14.1 is Associated with Pharmacokinetics and Pharmacodynamics of Bisoprolol.

Purpose: Bisoprolol is a widely used beta-blocker in patients with cardiovascular diseases. As with other beta-blockers, there is variability in response to bisoprolol, but the underlying reasons for this have not been clearly elucidated. Our aim was to investigate genetic factors that affect bisoprolol pharmacokinetics (PK) and pharmacodynamics (PD), and potentially the clinical outcomes. Patients and
Methods: Patients with non-ST elevation acute coronary syndrome were recruited prospectively on admission to hospital and followed up for up to 2 years. Patients from this cohort who were on treatment with bisoprolol, at any dose, had bisoprolol adherence data and a plasma sample, one month after discharge from index hospitalisation were included in the study. Individual bisoprolol clearance values were estimated using population pharmacokinetic modeling. Genome-wide association analysis after genotyping was undertaken using an Illumina HumanOmniExpressExome-8 v1.0 BeadChip array, while CYP2D6 copy number variations were determined by PCR techniques and phenotypes for CYP2D6 and CYP3A were inferred from the genotype. GWAS significant SNPs were analysed for heart rate response to bisoprolol in an independent cohort of hypertensive subjects.
Results: Six hundred twenty-two patients on bisoprolol underwent both PK and genome wide analysis. The mean (IQR) of the estimated clearance in this population was 13.6 (10.0-18.0) L/h. Bisoprolol clearance was associated with rs11029955 (p=7.17×10-9) mapped to the region of coiled-coil domain containing 34 region (CCDC34) on chromosome 11, and with rs116702638 (p=2.54×10-8). Each copy of the minor allele of rs11029955 was associated with 2.2 L/h increase in clearance. In an independent cohort of hypertensive subjects, rs11029955 was associated with 24-hour heart rate response to 4-week treatment with bisoprolol (p= 9.3×10-5), but not with rs116702638.
Conclusion: A novel locus on the chromosomal region 11p14.1 was associated with bisoprolol clearance in a real-world cohort of patients and was validated in independent cohort with a pharmacodynamic association.

Introduction

Beta-blockers are a widely prescribed class of drugs used in the management of several conditions including hypertension, angina pectoris, arrhythmias, heart failure with reduced ejection fraction, migraine prophylaxis, and as an adjunct in hyperthyroidism. Beta-blockers are also used in the secondary prevention of cardiovascular events following an acute coronary syndrome (ACS) where they have been shown to reduce mortality.1,2 The cardioprotective properties of beta-blockers are largely mediated by inhibiting the stimulatory effects of catecholamines on beta (β)1-adrenergic receptors which leads to a decrease in myocardial oxygen demand by reducing heart rate and stroke volume, suppression of the risk of ventricular arrhythmias, prolongs diastolic coronary perfusion, and reduces left ventricular remodelling.

Bisoprolol is a second generation, highly selective β1-adrenergic receptor antagonist with low affinity for β2-adrenergic receptors. Approximately 50% of bisoprolol is metabolised by cytochrome P450 (CYP) enzymes 2D6 and 3A4, with the remainder eliminated unchanged by the kidneys.3 The balanced elimination makes bisoprolol a safer beta-blocker for patients with impaired hepatic or kidney function, with a low risk of accumulation and consequently, a reduced risk of adverse reactions.

Genomic variation in CYP2D6 has been shown to affect the pharmacokinetic and pharmacodynamics of a number of drugs, including the β1-specific antagonist metoprolol.4 Indeed, implementation studies are underway to pre-emptively genotype CYP2D6 to guide metoprolol prescribing.5 By contrast, relatively little is known about the genomic factors that influence the pharmacokinetics and pharmacodynamics of bisoprolol.

The aims of this present study were two-fold: (a) to develop a Population (Pop) PK model for bisoprolol in patients following a non-ST elevation-acute coronary syndrome (NSTE-ACS) to identify any genetic factors determining clearance; and (b) determine whether the identified loci had any effects on heart rate, a PD parameter of beta-blockade. Two separate cohorts of patients were used to achieve these aims: PK were determined in a cohort of patients with NSTE-ACS in the Pharmacogenetics of Acute Coronary Syndrome (PhACS) study,6,7 while the PD parameter was determined in an independent cohort from the Genetics of Drug Responsiveness in Essential Hypertension (GENRES) study.8

Materials and Methods

Participants

A total of 1,470 patients hospitalised with NSTE-ACS were prospectively recruited into the PhACS study from 16 sites in the UK. They were prescribed bisoprolol on admission, at a dose determined by the clinician, and followed up for at least 12 months to assess clinical outcomes including further cardiovascular events and mortality, as described previously.6 The protocol was approved by Liverpool (Adult) Research Ethics Committee, UK, and informed consent was obtained from all study subjects in accordance with the Declaration of Helsinki.

Patients were included if they: i) were on treatment with bisoprolol, at any dose, at visit 2 (one month after discharge from index hospitalisation); ii) had bisoprolol adherence data recorded at visit 2; and iii) had a plasma sample collected at visit 2. Exclusion criteria were: i) bisoprolol plasma levels below the lower limit of quantification (LLOQ= 0.5 ng/mL); ii) the participant was excluded during pre-GWAS standard quality control (QC) procedures; iii) poor adherence to bisoprolol, defined as at least one pill missed in the week before visit 2 as assessed by the Brief Medication Questionnaire (BMQ);9 iv) dose, weight or time of blood collection was not recorded; and v) if they were categorised as an outlier in relation to bisoprolol plasma concentrations according to Tukey’s method within the R package jmuOutlier ().10

Quantification of Bisoprolol in Plasma at Steady State

Blood samples were collected in EDTA-containing vacutainer tubes at visit 2, and plasma was stored at −80°C. Bisoprolol concentrations were determined using a multiple-analyte HPLC-MS/MS assay developed and validated as described previously.7 Briefly, 50 µL of plasma spiked with deuterated internal standard (D5-bisoprolol) was extracted and analysed using the Shimadzu Nexera X2 modular system (Kyoto, Japan) coupled to a Sciex triple quadrupole 6500 QTRAP mass spectrometer (AB Sciex, Warrington, UK). Chromatographic separation was performed through a 2.7µm Halo C18 column with a gradient mixture of acetonitrile and water (10–90%) with 0.1% formic acid at a flow-rate of 500 µL/min to give a total run time of 6 minutes. The MS analysis was carried out in positive ionisation mode with multiple reaction monitoring (m/z 326.2 → 116.21 for bisoprolol and 331.0 → 121.1 for D5-bisoprolol). The assay was linear for the calibration range (0.5 to 125 ng/mL). The mean extraction recovery (CV) was 98.5 (11.0%) and 97.4 (7.0%) for the low and high QCs, respectively.

Population Pharmacokinetics Modelling

Sparse bisoprolol plasma levels at steady-state were utilised for nonlinear mixed-effects modelling using the Stochastic Approximation Expectation Minimization (SAEM) algorithm in Monolix 2018R1 (Lixoft, France). The base PK model was built as a one-compartment oral dose model with a fixed absorption rate constant (Ka = 2.31 h−1), fixed volume of distribution (V = 218 L), and the initial estimate for clearance (Cl) was 10.2 L/h. The model and the initial values were obtained from a previous bisoprolol PK study in patients with chronic heart failure.11

The interindividual variability of the PK parameter Cl was estimated using log-normal distributions with mean zero and variance ω2. In order to capture noise related to the analytical method, residual errors were calculated based on the inter-run standard deviation (SD) of the three independent runs. SD was inputted into the makeErrorPoly() function in Pmetrics package.12 This meant we could specify and fix the combined proportional and additive residual error model with means of zero and variances a and b. SDs (0.024, 0.068, 1.147 and 4.96) were obtained from six replicates at each of three independent runs for the following concentrations: 0.5, 1.5, 45 and 95 ng/mL.

In order to adjust for uncertainties in time after the last dose, the timing of every patient dose was set to 5AM to represent the hypothetical earliest time point that a patient would administer their bisoprolol, and the lag-time parameter (tlag) was utilised. The influence of the following covariates on the PK parameter Cl was assessed: body weight, sex, age, smoking status, use of diuretics, and chronic kidney disease (defined as a serum creatinine > 150µmol/L). The goodness of fit for nested models was compared using the difference in the −2xlog-likelihood (−2LL) in a forward selection process using the critical value of 3.84.

Genome-Wide Association Analysis

DNA samples obtained from PhACS participants were genotyped using the lllumina HumanOmniExpressExome-8 v1.0 BeadChip array at Edinburgh Genomics (Roslin Institute, Scotland). Standard per-individual QC filters were applied to the whole PhACS cohort (n=1,470) to exclude participants with: a sample call rate <95%, discordant clinical-genetic-derived gender, aberrant heterozygosity, one of each pair found to be cryptically related (identity by descent >0.1875), and genetically non-European participants based on the first two principal components (PC) from a PC analysis model built using international HapMap 3 data. Per-marker QC filters were similarly applied to the whole PhACS cohort to exclude single nucleotide polymorphisms (SNPs) with: call rate <95%, minor allele frequency <5%, or deviation from Hardy-Weinberg equilibrium (test p value <0.0001). These QC procedures excluded 113 participants, reducing the overall PhACS cohort to n=1,357, and led to a reduction in variants from 951,117 to 598,054. The genotype scaffold was pre-phased using SHAPEIT v213 and imputed to the 1000 Genomes Phase I reference panel14 using IMPUTE2.15 Post-imputation QC was performed by removing variants with a minor allele frequency <1% and info score <0.4.

The estimated mean individual clearance for bisoprolol was used as the phenotype in the GWAS. The GWAS was undertaken using multivariable linear regression in SNPtest version 2.4.1 (University of Oxford, Oxford, UK) under a frequentist additive dosage model, and adjusted for the first two PCs. A genome-wide statistical significance threshold of 5.0×10−8 was applied. LocusZoom16 was used to visualise the regions with genome-wide significant p-values.

CYP3A and CYP2D6 Genotype to Phenotype Analysis

As bisoprolol metabolism is mediated by CYP2D6 and CYP3A4, we performed a separated analysis focused on these genes. The following pharmacogenetic variants known to affect metabolic function of the enzymes were extracted from our GWAS: 1)CYP2D6 variants, *3 (rs35742686), *4 (rs3892097), *9 (rs5030656), *10 (rs1065852) and *41 (rs28371725), and 2)CYP3A variants, CYP3A4*22 (rs35599367) and CYP3A5*3 (rs776746). For imputed GWAS dataset, the values were categorized as follows: wild-type (WT, <0.3), heterozygous (0.7–1.3), and homozygous (>1.7). Participants with imputed genotype values outside these ranges were not included in the subsequent analysis. The number of CYP2D6 gene copies were identified with commercially available TaqmanTM copy number assay (ThermoFisherScientific). qPCR data were acquired using Quantstudio 6 Flex Real-Time PCR System (ThermoFisherScientific) and analysed with CopyCaller v2.1 software (LifeTechnologies). CYP2D6 alleles were grouped according to their function activity scores (AS), and phenotypic status was assigned as poor metabolizer (PM, AS=0), intermediate metabolizers (IM, AS=0.25–1), extensive metabolizers (EM, AS=1.25–2.25), or ultra-rapid metabolizers (UM, AS>2.25), using the CYP2D6 genotype-phenotype translation17 developed by the Clinical Pharmacogenetics Implementation Consortium and Dutch Pharmacogenetics Working Group. Genotype-phenotype translation for CYP2D6 and CYP3A present among PhACS participants is represented in the . The T allele for CYP3A4*22 (rs35599367)18 and GG genotype for CYP3A5*3 (rs776746)19 result in decreased enzymatic activity. We defined the metabolic status (phenotype) for CYP3A as described previously:20 PM were defined as at least one decreased activity allele (GG for rs776746 and CT or TT for rs35599367 - CYP3A5*3/*3* CYP3A4*1/*22 or CYP3A4*22/*22); IM with no CYP3A5 activity and no decreased activity allele in CYP3A4 (GG for rs776746 and CC for rs35599367 – CYP3A5*3/*3* CYP3A4*1/*1); or at least one functional allele in CYP3A5 and at least one decreased activity allele in CYP3A4 (AG or AA for rs776746 – CYP3A5*1/*3 or CYP3A5*1/*1 - and CT or TT for rs35599367 - CYP3A4*1/*22 or CYP3A4*22/*22); EM with at least one functional allele in CYP3A5 (AG or AA for rs776746 - CYP3A5*1/*3 or CYP3A5*1/*1) and no decreased function allele in CYP3A4 (CC for rs3559367 - CYP3A4*1/*1).

Validation Study

We were unable to identify another study which had evaluated the pharmacokinetics of bisoprolol with accompanying genetic data or samples to serve as a replication cohort. We therefore collaborated with the investigators in the GENRES (Genetics of Drug Responsiveness in Essential Hypertension) study.8 GENRES was a prospective, placebo-controlled trial, where men (age, 35–60 years) with mild-to-moderate hypertension received 4-week monotherapies with four antihypertensive drugs separated by 4-week placebo periods. Accordingly, a total of 216 moderately hypertensive men received bisoprolol 5mg/day as a monotherapy for 4 weeks. Our aim was to determine if any SNPs with p-values <5x10−6 in our PhACS study showed an association with heart rate response to bisoprolol in the GENRES study.

The DNA samples were genotyped using the Illumina HumanOmniExpress BeadChip (Illumina, San Diego, CA, USA).8 Phasing and imputation of the genotypes was done utilizing a Finnish population-specific reference panel of 3,775 high-coverage whole-genome sequences. Association analyses were undertaken with office and ambulatory 24-hour heart rate (HR) change on bisoprolol using mean values from up to four placebo periods as the baseline level, with an additive genetic model in Plink2 ().21 The baseline HR level (from either office or ambulatory 24-hour measurements, as appropriate) was used as a covariate. Other tested variables (age, body mass index, creatinine clearance and smoking) were not significant in stepwise regression analysis and were not included as covariates.

Results

Study Population and Plasma Drug Measurements

Of the 1,470 PhACS participants, a total of 653 met the eligibility criteria and were included in this study (). The demographic and clinical profiles of this cohort are shown in Table 1. As expected, plasma bisoprolol steady-state concentrations varied according to the daily bisoprolol prescribed dose (from 1.25 to 20 mg/day) ().

Table 1

Demographic and Clinical Profile of the Bisoprolol Pharmacokinetics Cohort

Demographics
Age (years, Mean ±SD)	64±12
Sex (% Female)	25.1
Weight (Kg, Mean ±SD)	83±16
BMI (Kg/m2, Mean ±SD)	29±6
Medical History
Hypertension (%)	53.3
Hyperlipidemia (%)	50.5
Diabetes (%)	18.8
Chronic renal failure (%)	5.4
Current smokers (%)	25.9
NYHA classification
I (N)	381
II (N)	204
III or IV (N)	49
IV (N)	5
Medications
Aspirin (%)	94.0
P2Y12 inhibitors (%)	83.9
Statins (%)	96.0
ACE inhibitors (%)	74.1
ARB (%)	8.3
Aldosterone antagonists (%)	3.2
Warfarin (%)	4.4
Calcium channel antagonists (%)	10.3
Potassium channel activators (%)	10.4
Diuretics (%)	20.4

Abbreviations: NYHA, New York Heart Association; ACE, angiotensin converting enzyme; ARB, angiotensin II receptor blockers.

The concentration-time relationship for bisoprolol was best described as a one-compartment model with first-order absorption and including tlag. The combined residual error model was calculated and fixed at 0.24 and 0.05, respectively. The PK model selection process included clinical covariates associated with bisoprolol clearance (). The final model is shown in Table 2, and included age, body weight, sex, smoking, and concomitant use of a diuretic, selected due to a lower −2LL value. Chronic kidney disease did not improve the model, and therefore was not included in the final model. Individual predictions versus observations, and prediction distribution of the final model are shown in the . The mean [25–75% interquartile range] of the estimated clearance in the studied population was 13.6 [10.0–18.0] L/h. The standard deviation of the random effects (omega) for clearance in the base model was estimated to be 0.455 (coefficient of variation 48.0%), and it was reduced to 0.402 (coefficient of variation 42.7%) in the final model with the covariates, showing that approximately 5% of the variability in clearance is explained by the clinical covariates.

Table 2

Pharmacokinetic Population Parameters of the Final Model

Parameter	Estimate	SE	p-value
Tlag_pop	3.17	0.149	-
ka_pop	2.31	FIXED	-
V_pop	218.00	FIXED	-
Cl_pop	13.1	0.382	-
beta_Cl_tAGE	−0.68	0.12	1.32x 10−8
beta_Cl_tWT	0.68	0.11	8.96x 10−10
beta_Cl_DIURETIC	−0.157	0.047	0.000882
beta_Cl_SMOKING	0.173	0.048	0.000335
beta_Cl_SEX	−0.122	0.047	0.00971
omega_Tlag	0.309	0.038	-
omega_ka	2.38	0.23	-
omega_V	0.193	0.03	-
omega_Cl	0.402	0.015	-
a	0.24	FIXED	-
b	0.05	FIXED	-

Abbreviations: Tlag, lag time before absorption; ka, absorption constant rate; V, volume of distribution; Cl, clearance, WT, body weight.

Genome-Wide Association Study for Clearance of Bisoprolol

Thirty-one samples were excluded due to the genetic QC procedures, leaving 622 individuals in the genome-wide association analysis. The Manhattan plot for the GWAS of bisoprolol clearance (Figure 1) revealed two genome-wide significant signals. The distributions of polymorphisms were in Hardy–Weinberg equilibrium. No genomic inflation was observed in the GWAS (λ= 0.999). The strongest association was observed in the region of coiled-coil domain containing 34 region (CCDC34) on chromosome 11. The lead SNP was rs11029955 (p=7.17x10−9) with five SNPs being in high linkage disequilibrium (r2 >0.8) with rs11029955. Each copy of the minor allele T (MAF=0.32) of rs11029955 was associated with an increase in clearance of bisoprolol by 2.2 L/h. The second signal, rs116702638 (p=2.54x10−8), mapped to a locus near LON peptidase N-terminal domain and ring finger 2 (LONRF2) on chromosome 2 (Table 3). Each copy of its minor allele (MAF=0.01) was associated with a 12.5 L/h increase in bisoprolol clearance.

Figure 1

Manhattan plot of bisoprolol clearance in patients with acute coronary syndrome. The red line represents the genome-wide significance p value threshold (5x10−8) and blue lines represents the nominal significance threshold (1x10−5).

Table 3

SNPs That Reached Genome-Wide Significance (p value < 5×10−8) for Bisoprolol Clearance in PhACS Patients

Variation	Chr	Position (Build 37)	Ref. Allele	Variant Allele	Variant Allele Frequency	Beta	P value	Info Score
rs11029955	11	27,354,487	G	T	0.32	2.22	7.17E-09	0.997
rs7935021	11	27,373,436	C	T	0.31	−2.18	1.12E-08	1.000
rs61887779	11	27,364,260	G	A	0.28	2.28	1.17E-08	0.992
rs1871254	11	27,369,906	T	C	0.28	2.27	1.20E-08	0.996
rs2291022	11	27,370,588	T	C	0.28	2.27	1.21E-08	0.996
rs116702638	2	100,915,330	G	T	0.01	12.54	2.54E-08	0.568

Note: Positive beta values indicate increase in clearance rate associated with the variant allele.

Abbreviations: Chr, chomosome; Ref, reference.

CYP3A and CYP2D6 Genotypes and Phenotypes and Clearance of Bisoprolol

No significant associations in the CYP3A and CYP2D6 regions were observed using the GWAS data. Regional plots for these loci are shown in Figure 2. Genotype-based phenotypes were assigned for CYP2D6 and CYP3A according to the presence of variant alleles and copy number variants (for CYP2D6). The genotype and phenotype frequencies are presented in the , respectively. We found no significant association of single CYP2D6 and CYP3A genotypes nor metabolizer phenotypes with clearance of bisoprolol. Furthermore, we found no significant differences between PM (12.88 [8.75–16.63] L/h, median [IQR], n=4) and EM (14.50[10.74–19.00] L/h, median [IQR], n=27) for both CYP3A and CYP3D6 enzymes.

Figure 2

Regional association plot for bisoprolol clearance at the CCDC34 locus (A), LONRF2 (B), CYP2D6 (C), and CYP3A (D). The colours reflect linkage disequilibrium (r2) for the signals with the lowest p values.

In the validation study (GENRES), we found a significant association of our lead SNP (rs11029955) with HR change on bisoprolol using ambulatory 24-hour HR measurements (p= 9.3x10−5, beta= −1.8) and office measurements (p= 0.049, beta=−1.1). However, no association was found with changes in office or ambulatory systolic and diastolic pressure after treatment with bisoprolol. No association was found with the SNP rs116702638 (Table 4). Results for variants with p values ≤ 5×10−6 in PhACS are shown in the .

Table 4

Associations Between Office and 24-Hour Ambulatory Changes in Heart Rate and the Identified Bisoprolol GWAS-Significant SNPs in GENRES Patients

			24h Ambulatory Measurements		Office Measurements
SNP	Variant Allele	VAF Allele Frequency	Beta	P value	Beta	P value
rs11029955	T	0.32	−1.8	9.3x10−5	−1.1	0.049
rs116702638	T	0.02	0.6	0.670	0.9	0.590

Notes: A negative beta value indicates a larger decrease in heart rate associated with the coded allele. The magnitude of beta indicates the magnitude of HR change per one copy of the coded allele. 206 individuals were analysed for 24h ambulatory measurements and 216 for office measurements.

Abbreviations: VAF, variant allele frequency. SNP: single nucleotide polymorphism.

Discussion

A recent systematic review and meta-analysis of pharmacogenetic studies undertaken with bisoprolol did not show convincing evidence of associations with either PK or PD genes, with the authors concluding that there is a “need of further studies about the impact of genetic variants on bisoprolol response”22. We have risen to this challenge, and to the best of our knowledge, our study represents the first GWAS of bisoprolol PK, based on bisoprolol clearance estimated by Pop-PK modeling. Importantly, two genome-wide significant loci were identified, with the lead SNP rs11029955 also shown to be associated with HR reduction in an independent validation cohort comprising individuals with hypertension. Interestingly, no associations were found using individual pharmacogenetic variants, nor genotype-based metabolizer phenotypes, in the two main enzymes involved in bisoprolol metabolism, CYP2D6 and CYP3A, consistent with the systematic review of Castano-Amores et al.22

Steady-state plasma concentrations of bisoprolol were obtained in 653 patients with a NSTE-ACS at a median of one month after hospital discharge, and PK parameters were modelled using nonlinear mixed-effects. The best Pop-PK model was a one-compartment oral dose model including tlag. The covariates age, body weight, use of a diuretic, smoking status, and sex, were significantly associated with bisoprolol clearance and were thus included in the final model. Those covariates explained approximately 5% of the variability in clearance in the model. Interestingly, a previous study has shown that smoking has a significant effect on bisoprolol clearance, possibly due to induction of CYP3A4.23 Women have higher drug exposure and lower clearance of metoprolol,24 but the influence of sex in the pharmacokinetics of bisoprolol has not been studied. The typical clearance estimated in our population was 13.1 L/h. Previous Pop-PK studies in healthy individuals, and patients with chronic heart failure, have reported mean bisoprolol clearances of 15.6 L/h3 and 10.2 L/h.11 The clearance observed in our study is in alignment with previous findings despite the differences in patient populations studied, study design, sampling times and modelling algorithms, all of which may influence the results. In our model, older age and lower body weight were most strongly associated with reduced bisoprolol clearance. It was previously reported that body weight and creatinine clearance were both associated with bisoprolol clearance in a middle-aged and elderly population.25 In our study, chronic kidney disease (a binary variable defined as a serum creatinine > 150µmol/L) did not improve the PK model. Actual creatinine clearance values were not available in PhACS participants and represent a limitation in the Pop-PK model described in this study. Population pharmacokinetic studies on bisoprolol are scarce and of limited interpretability due to differences in study designs and methodology. It has been suggested that large and well-designed population pharmacokinetics studies linked with genetic studies should be increasingly used in the future to identify factors associated with the pharmacokinetics variability of beta-blockers.26 This is exactly what we have done in this study, with the large sample size (n=653) being a strength.

Bisoprolol has a balanced elimination: 50% is enzymatically metabolised and 50% is eliminated unchanged in the urine. CYP2D6 and CYP3A4 are the main hepatic enzymes involved in the catalytic metabolism to O-deisopropylation of the drug. In vitro studies have demonstrated that CYP2D6 has a higher catalytic activity for bisoprolol compared to CYP3A4.27 However, due to the higher abundance of CYP3A4 in human microsomes, the contribution of CYP3A4 and CYP2D6 to the total metabolism of bisoprolol is estimated to be about 28% and 7%, respectively. The CYP3A subfamily comprises four isoenzymes (3A4, 3A5, 3A7, 3A43), and CYP3A4 and CYP3A5 have similar substrate specificities. Studies that have investigated the impact of genomic variation on bisoprolol PK are scarce. One small study with 40 patients with cardiovascular disorders showed that peak and trough bisoprolol plasma concentrations were similar in individuals with 0, 1, or 2 copies of the CYP2D6*10 reduction-of-function allele.28 According to the guidelines of the Dutch Pharmacogenetics Working Group, there are currently no dosing recommendation for bisoprolol based on CYP2D6 genotype (). In our study, bisoprolol clearance was not affected by pharmacogenetic variants in CYP3A4, CYP3A5 and CYP2D6 genes (including CYP2D6*10). CYP2D6 is highly polymorphic with more than 100 different alleles. We assigned CYP2D6 metabolizer status by combining CYP2D6 structural variants (presence of deletion or duplication) and five CYP2D6 variants most relevant to drug metabolism (*3, *4, *9, *10, and *41). We found that there was no significant association with any of the CYP2D6 variants we included in the study; additionally, there was no association with metabolizer phenotype. Similarly, we did not find an association with CYP3A4 and CYP3A5 allelic variants. Clearly it is possible that our findings are explained by the minor role of hepatic metabolism in bisoprolol clearance. Alternatively, bisoprolol exposure could be affected by genetic variants not captured in our GWAS platform, including rare variants, which is a limitation of our analysis.

We found two new genomic loci that showed a genome-wide association with estimated total clearance. The first is near the gene CCDC34 (coiled-coil domain containing 34 region) on chromosome 11. CCDC34 is a member of the coiled-coil domain containing and that is expressed in several human tissues. High expression of CCDC34 has been found in bladder cancer29 and colorectal cancer,30 and it has been suggested to play a role in cell proliferation and anti-apoptosis, and in promoting invasion and metastasis of cancer cells. The lead SNP at this locus, rs11029955, is in linkage disequilibrium with 4 other intronic CCDC34 SNPs (rs7935021, rs61887779, rs1871254, and rs2291022). Furthermore, rs1871254 has a rank of 5 on RegulomeDB, and so may affect binding of transcription factors,31 but was not an eQTL on analysis of the GTEx database. Unfortunately, very little information on signalling pathways is currently available, and the SNP is not associated with any disease or trait in the Genome-Wide Association Study catalogue.32 Interestingly, evaluation of the PheWAS database developed from >460,000 UK Biobank participants () showed that the G allele was associated with the medication code Cardicor® (bisoprolol) (beta= −0.00025, p=2.2 x10−4, based on 415 cases) and hypovolemic shock (beta = −0.00004, p=4.2 x 10−5, based on 10 cases). The SNP rs11029955 is present at relatively high frequencies in European (MAF 0.32), African (MAF 0.49), and Asian (0.42–0.50) populations, suggesting that the association may also be relevant for other ethnic groups. However, it is also important to state that we have not identified the causal variant, and it is possible that this is located in another gene in the vicinity, rather than in CCDC34. Functional studies will be needed to understand the mechanistic basis for our findings.

We could not find another cohort of patients with ACS treated with bisoprolol, where there was determination of both bisoprolol pharmacokinetics and genomic analysis, as a replication dataset. We therefore investigated the GENRES cohort and evaluated the effect of the lead SNPs on the effect of bisoprolol on heart rate, a pharmacodynamic parameter. Intriguingly, the T allele of our lead SNP (rs11029955), which was associated with increased bisoprolol clearance, was associated with a larger decrease in HR after four weeks on bisoprolol in hypertensive individuals in GENRES, but there was no association with blood pressure. It is well known that systemic drug exposure is an important determinant of drug efficacy, and thus variability in plasma exposure can lead to variability in PD response. Based on the observed increase in bisoprolol clearance associated with the rs11029955 variant allele, we would have expected a smaller reduction in HR in individuals in GENRES carrying the variant allele, but the opposite was observed. However, we must acknowledge that there are important differences between the 2 cohorts, and GENRES cannot be considered a true replication cohort for PhACS. Whilst patients in PhACS had multiple comorbidities and were on multiple drugs, GENRES participants were all males and had mild-to-moderate hypertension with no clinically significant complications. Moreover, variability in disease state, a drug’s molecular target and downstream signalling pathways can also lead to variable drug response independent of drug PK. We therefore speculate that the genomic region may regulate pathways involved in cardiovascular and renal homeostasis, and that increased clearance may be a protective mechanism to mitigate the PD effects of bisoprolol. This is also consistent with the PheWAS finding that the SNP was associated with the use of bisoprolol and hypovolemic shock. Further studies are required to understand the function of CCDC34 in the body.

The second signal identified in our GWAS, rs116702638, and mapped to LON peptidase N-terminal domain and ring finger 2 locus (LONRF2) on chromosome 2. LONRF2 is expressed at higher levels in adrenal, brain, and endometrium compared to other human tissues.33 However, this is likely to be a false positive as the frequency of the SNP is very low in the population, and we could not find any association in the validation cohort.

The systematic review of bisoprolol pharmacogenetics showed that β-adrenoceptor genotypes (ADRB1 and ADRB2) have been the most widely studied, but the results have been contradictory.22 No association was found with these genes in our study. Interestingly in the GENRES study,8 polymorphisms in ACY3 (aminoacylase III) on chromosome 11 were significantly associated with the extent of blood pressure lowering with bisoprolol. However, the lead SNP within ACY3 is 40 Mb away from our GWAS-significant SNP, rs11029955, and they are not in linkage disequilibrium. It is possible that there may be a causal locus on chromosome 11 between the locus identified in our study and the ACY3 locus, and this will require further in-depth investigation. Blood pressure response to beta-blockers (metoprolol, atenolol and bisoprolol) in patients of European Ancestry has also been linked to a missense variant in BST1 (chromosome 4) by the International Consortium of Antihypertensive Pharmacogenomics Studies,34 but whether this is drug-specific or class-specific will need further study.

The main strength of our study is its use of a large real-world patient cohort (with associated drug concentrations), which reflects the clinical practice setting of cardiovascular secondary prevention. One of the limitations, however, was that the time since last bisoprolol dose was estimated, because only the time of blood sampling was precisely recorded in PhACS. However, this represents a common limitation in real world observational pharmacogenomic studies,35 that cannot be retrospectively collected. Importantly, we have addressed this uncertainty using the tlag parameter in the Pop-PK modelling.

Conclusion

Our GWAS showed that the SNP rs11029955 was associated with both bisoprolol clearance in our PhACS cohort, and with HR reduction in the independent GENRES cohort. No associations were found with CYP2D6 or CYP3A genotypes or metabolizer status despite the fact that these two P450 enzymes are involved in bisoprolol metabolism. Further studies are required to understand the functional effects of rs11029955, its role in cardiovascular and renal homeostasis and therefore on bisoprolol PK/PD.