Analysis of microRNA processing machinery gene (DROSHA, DICER1, RAN, and XPO5) variants association with end-stage renal disease.

BACKGROUND: MicroRNA (miRNA) processing machinery gene variant was associated with several diseases. We aimed to explore for the first time the association of machinery gene (DROSHA rs10719A/G; DICER1 rs3742330A/G; RAN rs14035C/T; and XPO5 rs11077T/G) variants with the susceptibility and phenotype of end-stage renal disease (ESRD).
METHOD: A total of 281 participants (98 ESRD patients and 183 healthy volunteers) were enrolled. Real-Time TaqMan allelic discrimination assay was applied for the genotyping of the specified variants. Multiple logistic regression models, univariate, multivariate, and principal component analyses were carried out.
RESULTS: Carrying one DICER1 rs3742330*G allele conferred protection against developing ESRD [heterozygote comparison: OR = 0.30, 95% CI = 0.15-0.62, dominant model: OR = 0.35, 95% CI = 0.17-0.70]. Similarly, for XPO5 rs11077T/G, homozygote and heterozygote carriers of G variant were less likely to develop ESRD [homozygote comparison: adjusted OR = 0.23, 95% CI = 0.11-0.50, and heterozygote comparison: OR = 0.50, 95% CI = 0.22-0.92, and allelic model: OR = 0.46, 95% CI = 0.34-0.70]. RAN rs14035*TT subjects were 5 times more likely to develop ESRD while being heterozygote (CT) have twice the risk [homozygote comparison: 5.18, 95% CI = 2.28-11.8, heterozygote comparison: OR = 2.12, 95% CI = 1.10-409]. Subgroup analysis also detected DICER1 rs3742330*A, XPO5 rs11077*T, and RAN rs14035*T as genetic risk determinants for ESRD development in both sex and age categories. Two genotype combinations of DROSHA/DICER1/XPO5/RAN were associated with increased susceptibility to ESRD [A-A-T-T: OR = 9.49, 95%CI = 2.48-36.31 (P = .001) and G-A-T-T: OR = 5.92, 95%CI = 1.77-19.83 (P = .004), respectively].
CONCLUSION: Variants and combined genotypes of DICER1 rs3742330, XPO5 rs11077, and RAN rs14035 might be associated with ESRD development in the study population. Integrating molecular analysis in ESRD risk stratification is warranted.

INTRODUCTION

Despite advancements in health care, mortality and morbidity rates are still high among end‐stage renal disease (ESRD) cases. , Hence, searching and applications for early identification and treatment strategies for this devastating disease are in need worldwide. Currently, many genome‐wide association studies have suggested a potential association between single‐nucleotide polymorphisms (SNPs) and a quantitative trait or disease of interest. ,

Several studies have proved microRNAs’ (miRNAs) implication in renal development, homeostasis, and physiological functions. , , , These types of short non‐coding RNAs (21‐25 nucleotides) were reported to regulate virtually all major cellular processes by inhibiting target gene expression at the post‐transcriptional level through elicit either sequence‐specific translational repression or messenger RNA (mRNA) degradation. Consequently, miRNAs are considered influential in the outcome of various cellular activities under physiological conditions (including normal kidney function and homeostasis) and kidney diseases. , , miRNAs’ biosynthesis involves several miRNA machinery genes and occurs in multiple steps. “RNA polymerase II” produces primary miRNA (500‐3000 nucleotides), which is processed by a multiprotein complex that includes DROSHA, a Class 2 ribonuclease III enzyme, to form a precursor miRNA hairpin (60‐100 nucleotides; pre‐miRNA), After pre‐miRNA has been exported to the cytoplasm by Ran GTPase (RAN) and exportin 5 (XPO5), it is further processed by DICER1, a polymerase II enzyme. A double‐stranded miRNA duplex that unwinds and forms a single‐stranded mature miRNA is generated. Despite few miRNA SNPs have been implicated in chronic kidney disease risk, in part, through impact miRNA processing with a subsequent reduction in gene expression, , similar evidence for genes involved in miRNA biosynthesis machinery is lacking. In this sense, the authors were inspired to investigate for the first time the possible association between the miRNA machinery gene DROSHA rs10719; DICER1 rs3742330; RAN rs14035; and XPO5 rs11077 variants with the risk of ESRD in the present population, and to correlate the findings with the available clinical and laboratory data of the patients. The criteria and the selection process of these specified variants are detailed in the Section 2.3 of this study. The identification and characterization of such associations may highlight their roles in disease susceptibility and help in the development of novel genetic risk stratification for targeted screening and management in the near future.

SUBJECTS AND METHODS

Study participants

The present case‐control study enrolled 98 unrelated ESRD patients (36 males and 62 females) undergoing regular hemodialysis treatment (3 times a week) in the Nephrology Center of “Mohammed bin Saud Al‐Kabeer for renal dialysis,” Arar central hospital, Northern Border area, Saudi Arabia. Patients aged >18 years on regular hemodialysis and with no history of severe diseases as malignancy. Evidence of ESRD by consistently elevated serum creatinine level >4.0 mg/dL and/or estimated glomerular filtration rate (eGFR) <15 mL/min/1.73 m2 was confirmed from the patient records according to the local center protocols. For each patient, the basic information such as age, sex, blood pressure, urinary protein level, blood urea nitrogen, and complete lipid profile were retrieved. One hundred eighty‐three ethnically matched Saudi individuals from the same geographic area attending the “Out‐patient Care Unit” during their routine health check were considered as the control group. This group had no history or laboratory findings of chronic kidney diseases (ie, estimated glomerular filtration rate (eGFR) >60 mL/min/1.73 m2)). The study was conducted following the guidelines in the Declaration of Helsinki and was approved by the local Medical and Bioethics committee (approval No. MED‐2017‐1‐8‐F‐7381). Informed consent was obtained from each patient before taking part.

Sample collection and laboratory analysis

Five‐milliliter fasting blood was collected from all participants (for patients, the samples were withdrawn before the second dialysis session of the week) on plain vacutainer serum separator tubes and EDTA tubes. The separated serum from the former tubes after centrifugation was subjected to the routine hospital schedule of biochemical analyses, including the kidney function tests. The EDTA tubes were used for DNA extraction. The available results of the biochemical assays that done routinely for the patients were taken from their records. Otherwise, collected sera from the control group were tested for all routine measurements, as explained in our previous work. All the quality control measurements of the laboratory work were applied according to the supplier instructions and the local protocols.

Criteria for selecting the study variants and in silico data analysis

Literature review and in silico approach were applied for retrieving the microRNA biogenesis pathway and selection of common (ie, minor allele frequency >0.05) SNPs located in miRNA machinery genes coding for four proteins related to microRNA processing and maturation. , , , Genomic sequence and variants were analyzed in ensemble.org, and common variants were selected for each gene (Table 1). Subcellular localization was identified in the COMPARTMENTS database (compartments.jensenlab.org). Protein‐protein interaction and gene ontology were analyzed using String v10.5 (string‐db.org).

Table 1

Selected microRNA machinery gene variants

Locus	Position	Gene	SNP ID	Alleles	MAF	Type
5p13.3	5:31401340	DROSHA	rs10719	A/G	0.48 (A)	3'UTR
14q32.13	14:95087025	DICER1	rs3742330	A/G	0.14 (G)	3'UTR
6p21.1	6:43523209	XPO5	rs11077	T/G	0.40 (G)	3'UTR
12q24.33	12:130876696	RAN	rs14035	C/T	0.27 (T)	3'UTR

Abbreviations: MAF, minor allele frequency; SNP, single‐nucleotide polymorphism. Data source: www.ensembl.org

MiRNAs machinery gene variants allelic discrimination analysis

Genomic DNA was extracted from whole blood using the QIAamp DNA Blood Mini kit (Qiagen) following the manufacturer's protocol. Extracted DNA purity and concentration were assessed by NanoDrop ND‐1000 (NanoDrop Technologies, Inc). Genotyping for the selected four SNPs in miRNA processing genes (DROSHA rs10719; DICER1 rs3742330; RAN rs14035; and XPO5 rs11077) were assayed using real‐time polymerase chain reaction (PCR) allelic discrimination technology as described in detail in our previous work. PCR reactions were run blindly in duplicates in a 25‐µL final volume containing 20 ng genomic DNA, TaqMan Universal PCR Master Mix, no UNG, and TaqMan SNP Genotyping Assay Mix according to the standard protocols. Appropriate controls were used in each reaction. PCR amplification was done using the StepOne Real‐Time PCR System (Applied Biosystems). Allelic discrimination was called by the SDS software version 1.3.1 (Applied Biosystems). Ten percent of the samples were re‐evaluated in separate runs with 100% concordance rate.

Statistical analysis

R 3.5.1, SPSS v23.0, and GraphPad Prism v7.0 were used for statistical analysis. Continuous variables were presented as mean and standard deviation (SD), while categorical data were expressed as frequency and percentage. Student's t, one‐way ANOVA, Mann‐Whitney U, Kruskal‐Wallis, and chi‐square tests were applied as appropriate. The P‐value was set to be significant at <.05. Genotype and allele frequencies and carriage rates between patients and controls were estimated as described previously. Hardy‐Weinberg equilibrium was calculated online (http://www.oege.org/software/hwe‐mr‐calc.shtml) and analyzed by the chi‐square test. Multiple SNP analysis was carried out under unconditional logistic regression models using the online tool SNPStats software (https://www.snpstats.net/start.html). Association between the four variants and ESRD risk was expressed as odds ratio (OR) with 95% confidence intervals (CI), adjusted for the covariates like age, sex, BMI, hypertension, diabetes, and the significant laboratory findings. Multiple logistic regression models (allelic, codominant, dominant, and recessive) were applied for genetic associations. Pearson's and Spearman's rank correlation analyses were employed. Linear regression models were run to assess the proportion of variation of albumin/creatinine ratios in the study population. The evaluation for the models was performed using collinearity diagnostics (tolerance and variance inflation factor). The principal component analysis was carried out to represent multivariate similarities and clustering using the “factoextra” and “FactoMineR” packages in R version 3.5.3 and R studio version 1.1.383.

RESULTS

Functional enrichment analysis

Drosha, ribonuclease III enzyme, is involved in the initial step of miRNA biogenesis. It cleaves primary miRNAs to form hairpin‐shaped precursor miRNAs in the nucleus. Exportin‐5 and the GTPase RAN proteins mediate their nuclear export to the cytoplasm where DICER1, another ribonuclease III, further trims pre‐miRNAs in the cytoplasm to generate miRNA duplex. Within the RNA‐induced silencing complex (RISC), mature miRNAs serve as a guide to complementary RNAs to degrade them or prevent their translation. Their subcellular localization and protein‐protein network are demonstrated in Figure 1.

Figure 1

In silico data analysis. (A‐D) subcellular localization of DROSHA, DICER1, exportin‐5, and RAN proteins. (E) Protein‐protein interaction (https://string‐db.org/cgi/network.pl?taskId=tkj6nrf3Olol). Abbreviations: DICER1: Double‐stranded RNA‐specific Endoribonuclease; DGCR8, DiGeorge Syndrome Critical Region Gene 8; DROSHA, Drosha Ribonuclease III; EIF2C1, Eukaryotic translation initiation factor 2C, 1; EIF2C2, Eukaryotic translation initiation factor 2C, 2; EIF2C4, Eukaryotic translation initiation factor 2C, 4; KPNB1, Karyopherin (importin) beta 1; NUTF2, Nuclear transport factor 2; RAN, RAN, Member RAS Oncogene Family; RANBP1, RAN‐binding protein 1; RANBP2, RAN‐binding protein 2; RCC1, Regulator of chromosome condensation 1; TARBP2, TAR (HIV‐1) RNA‐binding protein 2

Baseline characteristics of the study population

The demographic and clinical features of both ESRD patients and controls are depicted in Table 2. The mean age of patients was 46.6 ± 10.6 years. Nearly two‐thirds were females. Their BMI was significantly lower than controls. Diabetes (32.7%) and hypertension (57.1%) were common findings among patients. Laboratory results are presented in Figure 2. Statistically significant differences were found for almost all parameters.

Table 2

Demographic and clinical characteristics of the study population

		Controls (n = 183)	Cases (n = 98)	P‐value
Age, y	mean ± SD	43.7 ± 12.0	46.6 ± 10.6	.043
Age categories	<50 y	113 (617)	60 (61.2)	.931
	≥50 y	70 (38.3)	38 (38.8)
Sex	Female	86 (47.0)	62 (63.3)	.012
	Male	97 (53.0)	36 (36.7)
BMI, Kg/m2	mean ± SD	32.0 ± 7.0	30.1 ± 8.3	.039
Obesity	Negative	75 (41.0)	50 (51.0)	.131
	Positive	108 (59.0)	48 (49.0)
Hypertension	Negative	180 (98.4)	42 (42.9)	<.001
	Positive	3 (1.6)	56 (57.1)
Diabetes	Negative	182 (99.5)	66 (67.3)	<.001
	Positive	1 (0.5)	32 (32.7)

Data are presented as number (percentage) or mean ± standard deviation (SD).

Chi‐square and Student's t tests were used.

Statistically significant values at P < .05 are bold.

Abbreviation: BMI, body mass index.

Figure 2

Biochemical findings in ESRD patients and controls. Mann‐Whitney U test was used for comparison. P < .05 was considered statistically significant. Abbreviations: Alb/Cre, albumin creatinine ratio in urine; BUN, blood urea nitrogen; CRE, creatinine; Hb, hemoglobin; HDL, high‐density lipoprotein cholesterol; HS‐CRP, high specific C‐reactive protein; TC, total cholesterol; TG, triglyceride; URCA, uric acid

Single variant analysis of machinery genes

Apart from the XPO5 rs11077 gene variant, genotype frequencies of all SNPs showed deviation from that expected in HWE (Table S1). Minor allele frequencies in the current study compared to other populations worldwide are demonstrated in Tables S2. DROSHA rs10719*G and XPO5 rs11077*T had similar frequencies to African values, and DICER1 rs3742330*G showed consistent frequency with East Asians, while RAN rs14035*T value was close to that of Africans, Americans, and Europeans.

Genotype and allele frequencies and carriage rates are shown in Table 3. Multiple logistic regressions were applied to determine the influence of each SNP on disease risk under different genetic inheritance models (Table 4).

Table 3

Genotype and allele frequencies of variants in ESRD cases and controls

		Controls (n = 183)	Cases (n = 98)	P‐value	Crude OR (95%CI)	Adjusted OR (95%CI)
DROSHA rs10719 (A/G)
Genotype frequency	AA	62 (33.9)	23 (23.5)	.134	Reference	Reference
	AG	72 (39.3)	40 (40.8)		1.49 (0.80‐2.77)	1.39(0.74‐2.62)
	GG	49 (26.8)	35 (35.7)		1.92 (1.00‐3.67)	1.84(0.95‐3.58)
Allele frequency	A	196 (53.6)	86 (43.9)	.028	Reference
	G	170 (46.4)	110 (56.1)		1.47 (1.04‐2.09)
Carriage rate	A	134 (73.2)	63 (64.3)	.192	Reference
	G	121 (66.1)	75 (76.5)		1.31 (0.86‐2.00)
DICER1 rs3742330 (A/G)
Genotype frequency	AA	18 (9.8)	23 (23.5)	.004	Reference	Reference
	AG	162 (88.5)	75 (76.5)		0.36 (0.18‐0.71)	0.30(0.15‐0.62)
	GG	3 (1.6)	0 (0.0)		0.11 (0.01‐2.31)	0.00(0.00‐NA)
Allele frequency	A	198 (54.1)	121 (61.7)	.081	Reference
	G	168 (45.9)	75 (38.3)		0.73 (0.51‐1.04)
Carriage rate	A	180 (98.4)	98 (100)	.335	Reference
	G	165 (90.2)	75 (76.5)		0.83 (0.57‐1.20)
XPO5 rs11077 (T/G)
Genotype frequency	TT	27 (14.8)	30 (30.6)	.001	Reference	Reference
	TG	77 (42.1)	44 (44.9)		0.51 (0.27‐0.97)	0.50 (0.22‐0.92)
	GG	79 (43.2)	24 (24.5)		0.27 (0.13‐0.54)	0.23 (0.11‐0.50)
Allele frequency	T	131 (35.8)	104 (53.1)	<.001	Reference
	G	235 (64.2)	92 (46.9)		0.49 (0.34‐0.70)
Carriage rate	T	104 (56.8)	74 (75.5)	.019	Reference
	G	156 (85.2)	68 (69.4)		0.61 (0.40‐0.92)
RAN rs14035 (C/T)
Genotype frequency	CC	61 (33.3)	16 (16.3)	<.001	Reference	Reference
	CT	101 (55.2)	56 (57.1)		2.11 (1.11‐4.00)	2.12(1.10‐4.09)
	TT	21 (11.5)	26 (26.5)		4.72 (2.12‐10.4)	5.18(2.28‐11.8)
Allele frequency	C	223 (60.9)	88 (44.9)	<.001	Reference
	T	143 (39.1)	108 (55.1)		1.91 (1.34‐2.72)
Carriage rate	C	162 (88.5)	72 (73.5)	.039	Reference
	T	122 (66.7)	82 (83.7)		1.51 (1.01‐2.24)

Values are presented as number (percentage).

ESRD: end‐stage renal disease. Chi‐square test was used.

Statistical significance at P < .05 are in bold.

Age, sex, body mass index, hypertension, diabetes, and significant laboratory findings were used for adjusted OR calculations.

Table 4

ESRD risk under different genetic inheritance models

SNP	Allelic model	Homozygote comparison	Heterozygote comparison	Dominant model	Recessive model
rs10719	1.47 (1.04‐2.09)	1.91 (1.00‐3.69)	1.49 (0.80‐2.79)	1.66 (0.95‐2.95)	1.51 (0.89‐2.57)
rs3742330	0.73 (0.51‐1.04)	0.11 (0.01‐2.31)	0.36 (0.18‐0.71)	0.35 (0.17‐0.70)	0.26 (0.01‐5.12)
rs11077	0.49 (0.34‐0.70)	0.27 (0.13‐0.55)	0.51 (0.27‐0.97)	0.39 (0.21‐0.71)	0.42 (0.24‐0.73)
rs14035	1.91 (1.34‐2.72)	4.65 (2.10‐10.5)	2.10 (1.12‐4.08)	2.55 (1.39‐4.85)	2.75 (1.46‐5.31)

Two‐sided chi‐square (χ2) test was used.

ESRD: end‐stage renal disease; OR (95% CI), odds ratio, and confidence interval.

Age, sex, body mass index, hypertension, diabetes, and significant laboratory findings were used for adjusted OR calculations.

Bold values are statistically significant at P < .05.

DROSHA rs10719 (A/G) analysis revealed that the rs10719*G allele was more prevalent in ESRD patients than in the control group (56.1% vs 46.4%, P‐value = .028). However, no significant difference was observed in genotype frequencies between the two groups (P = .13).

Concerning DICER1 rs3742330 (A/G) polymorphism, a significantly higher frequency of rs742330*AA genotype was detected in patients accounting for 23.5%, while the same genotype was represented in only 9.8% of controls (P = .004). On the other hand, carrying one G allele conferred protection against developing ESRD [heterozygote comparison (AG vs AA): OR = 0.30, 95% CI = 0.15‐0.62, dominant model (GG + AG vs AA): OR = 0.35, 95% CI = 0.17‐0.70].

For the XPO5 rs11077 (T/G) variant, the minor allele frequency of rs11077 (T) was the risk allele. It was significantly higher in ESRD subjects (0.53) as compared to controls (0.36) (P < .001). As compared with the control group, rs11077*GG genotype was less prevalent in patients (24.5%) than controls (43.2%) (P = .001). Homozygote and heterozygote carriers of G variant were less likely to develop ESRD [Homozygote comparison (GG vs TT): adjusted OR = 0.23, 95% CI = 0.11‐0.50, and heterozygote comparison (TG vs TT): adjusted OR = 0.50, 95% CI = 0.22‐0.92, and allelic model (G vs T): OR = 0.46, 95% CI = 0.34‐0.70].

By genotype and allele frequency analysis of RAN rs14035 (C/T), a significant difference between patient and control groups was observed. Increased frequency of rs14035*T allele was found in patients (55.1% vs 39.1%, P < .001) [allelic model (T vs C): OR = 1.91, 95% CI = 1.34‐2.72]. Similarly, RAN rs14035*TT genotype was more predominant in ESRD patients (26.5%) compared to controls (11.5%) (P < .001). TT subjects are five times more likely to develop ESRD while being heterozygote (CT) have twice the risk [Homozygote comparison (TT vs CC): 5.18, 95% CI = 2.28‐11.8, heterozygote comparison (CT vs CC): OR = 2.12, 95% CI = 1.10‐409].

Stratification analysis by age and gender

Like the overall analysis, stratification by sex and age revealed no significant difference in genotype frequencies of DROSHA rs10719 variant (A/G) between patients and controls in either sex or both age categories. In contrast, DICER1 rs3742330*A, XPO5 rs11077*T, and RAN rs14035*T represented genetic risk determinants for ESRD development in both sex and age categories (Figure S1). The interaction analysis of gene variants with covariates (age and six) is illustrated in Table S3.

Multiple variant analyses of machinery genes

Correlation analysis was performed to explore gene‐gene interactions (Figure S2). Positive weak correlation was found between DICER1 rs3742330 and XPO5 rs11077 variants (r = .16, P = .006). RAN rs14035 SNP was negatively associated with DROSHA rs10719 (r = −.12, P = .035), DICER1 rs3742330 (r = −.31, P < .001), and XPO5 rs11077 (r = −.21, P < .001).

Association analysis of combined genotypes with disease risk was performed using the web‐based SNPstats software (Table 5). Three genotype combinations (A‐G‐G‐C, A‐A‐T‐T, and G‐A‐G‐C) were most common in the study population. Two genotype combinations of DROSHA (A/G), DICER1 (A/G), XPO5 (T/G), and RAN (C/T) were associated with an increased susceptibility to have ESRD [A‐A‐T‐T: OR = 9.49, 95% CI = 2.48‐36.31 (P = .001) and G‐A‐T‐T: OR = 5.92, 95% CI = 1.77 ‐ 19.83 (P = .004), respectively].

Table 5

Genotype combinations and the risk of developing end‐stage renal disease

	DROSHA (A/G)	DICER1 (A/G)	XPO5 (T/G)	RAN (C/T)	Frequency	Adjusted OR (95%CI)	P‐value
1	A	G	G	C	0.147	1.00	‐‐‐
2	A	A	T	T	0.144	9.49 (2.48‐36.31)	.001
3	G	A	G	C	0.103	2.13 (0.61‐7.48)	.24
4	G	G	G	C	0.094	0.90 (0.19‐4.35)	.9
5	G	A	T	T	0.078	5.92 (1.77‐19.83)	.004
6	G	A	G	T	0.076	2.65 (0.73‐9.70)	.14
7	A	A	G	T	0.066	3.01 (0.53‐16.89)	.21
8	G	G	T	C	0.064	2.37 (0.66‐8.55)	.19
9	A	G	T	C	0.046	0.28 (0.04‐1.90)	.19
10	G	A	T	C	0.038	4.62 (1.00‐21.41)	.052
11	G	G	G	T	0.035	7.49 (1.00‐51.22)	.051
12	A	A	T	C	0.031	1.98 (0.21‐18.45)	.55
13	A	G	G	T	0.029	0.00 (−Inf‐Inf)	1
14	A	A	G	C	0.027	0.00 (−Inf‐Inf)	1
15	A	G	T	T	0.008	0.00 (−Inf‐Inf)	1

The log‐additive model for interaction between genotype combinations and covariates was performed.

Adjusted by age, sex, body mass index, hypertension, diabetes, and significant laboratory findings. Statistically significant values at P < .05 are bold.

Global association P‐value: <.001.

GGTT combination was absent in the study population

Univariate analysis

Association analysis of demographic and genetic variables is shown in Table 6. DICER1 rs3742330 SNP was not associated with any clinical or biochemical variables. Diabetes mellitus was highly frequent in XPO5 rs11077*T carriers (TT and TG) (P = .040). XPO5 heterozygosity rs11077*TG was more prevalent in men (52.6%) compared to homozygote genotypes (TT: 14.3% and GG: 33.1%) (P = .005). For laboratory results, the minor allele of DROSHA rs10719*G exhibited marked serum creatinine levels (P = .040) and lower urinary albumin/creatinine ratio (P = .016).

Table 6

Association and correlation of gene variants with the clinical and laboratory findings in the ESRD cohort

	DROSHA (A/G)		DICER1 (A/G)		XPO5 (T/G)		RAN (C/T)
	P Ass	r (P cor)	P Ass	r (P cor)	P Ass	r (P cor)	P Ass	r (P cor)
Age	.331	−.058 (.333)	.582	.054 (.366)	.695	−.029 (.633)	.555	.010 (.865)
Sex	.454	−.069 (.248)	.229	−.058 (.335)	.005	.031 (.602)	.440	−.020 (.744)
BMI	.062	−.137 (.202)	.643	−.049 (.412)	.449	−.078 (.192)	.360	.038 (.527)
Obesity	.124	−.106 (.077)	.319	−.075 (.209)	.297	−.093 (.121)	.409	.051 (.396)
HTN	.496	.070 (.242)	.578	−.047 (.428)	.224	−.103 (.086)	.565	.044 (.463)
Diabetes	.425	.044 (.458)	.685	−.046 (.455)	.040	−.141 (.018)	.110	.043 (.478)
Hb	.268	.102 (.170)	.432	−.108 (.145)	.499	.095 (.200)	.058	.069 (.348)
TC	.542	−.028 (.703)	.310	−.100 (.175)	.648	−.122 (.129)	.724	−.008 (.913)
TG	.062	−.010 (.898)	.901	−.007 (.929)	.536	−.013 (.857)	.087	.139 (.060)
LDL	.720	−.132 (.074)	.499	.032 (.665)	.987	.041 (.584)	.542	−.069 (.355)
HDL	.947	.017 (.824)	.913	−.004 (.961)	.362	−.094 (.205)	.402	−.035 (.635)
BUN	.784	.055 (0.455)	.887	−.027 (.712)	.780	.030 (.687)	.742	.055 (.456)
CRE	.040	.184 (.012)	.991	.013 (.859)	.063	.074 (.321)	.074	.034 (.646)
URCA	.405	.098 (.184)	.583	.079 (.288)	.382	.096 (.194)	.415	−.114 (.122)
HS‐CRP	.055	.107 (.144)	.628	.075 (.313)	.150	.150 (.401)	.581	−.123 (.096)
Alb/Cre	.016	−.211 (.004)	.679	−.058 (.437)	.532	−.009 (.904)	.158	−.017 (.823)

Data are presented as P‐values for association (PAss) or correlation coefficient (r) and its P‐value (Pcor). Chi‐square, Kruskal‐Wallis, and Pearson's correlation tests were applied. Statistically significant values at P < .05 are bold.

Abbreviations: Alb/Cre, albumin creatinine ratio in urine; BMI, body mass index; BUN: blood urea nitrogen, BUN, blood urea nitrogen; CRE, creatinine; Hb, hemoglobin; HDL, high‐density lipoprotein cholesterol; HTN, hypertension; HS‐CRP, high sensitivity C‐reactive protein; TC, total cholesterol; TG, triglyceride; URCA, uric acid.

Multivariate analysis

Multiple regression analysis was applied to identify the role of the four SNPs in albumin/creatinine levels. The dominant model for the interaction between genotype combinations and covariates is shown in Table 7. Differential impact of various combinations of machinery SNP alleles was observed.

Table 7

Regression analysis of the association between SNPs and urinary albumin to creatinine ratio.

	DROSHA (A/G)	DICER1 (A/G)	XPO5 (T/G)	RAN (C/T)	Frequency	Mean difference (95%CI)	P‐value
1	A	G	G	C	0.128	0.00	‐‐‐
2	A	A	T	T	0.124	NA (NA‐NA)	NA
3	G	G	G	C	0.106	‐60.11 (‐109.72‐−10.49)	.018
4	G	A	G	C	0.096	59.45 (0.29‐118.6)	.049
5	G	A	T	T	0.095	92.92 (51.41‐134.43)	<.001
6	A	A	G	T	0.076	−23.07 (−77.05‐30.9)	.400
7	G	A	G	T	0.073	−34.03 (−91.04‐22.98)	.240
8	A	G	T	C	0.064	−73.93 (‐114.28‐−33.57)	.0003
9	G	G	T	C	0.056	−75.44 (−105.59‐−45.28)	<.001
10	G	A	T	C	0.039	4.15 (−10.72‐19.03)	.580
11	A	A	G	C	0.038	−56.54 (−77.85‐−35.23)	<.001
12	A	G	G	T	0.031	31.36 (20.33‐42.4)	<.001
13	G	G	G	T	0.030	−42.72 (−54.63‐−30.8)	<.001
14	A	A	T	C	0.023	−40.38 (−55.72‐−25.03)	<.001
15	A	G	T	T	0.014	−70.53 (−73.87‐−67.18)	<.001

Bold values indicate significance at P <.05.

The principal component analysis was applied to visualize the variations and similarities between study subjects based on their clinical, biochemical, and molecular features (Figure 3). Exploration across different types of variables demonstrated the influence of hypertension, diabetes, high age, female six, and shortness in ESRD patients (Figure 3A). On the other hand, laboratory findings resulted in a better demarcation between patients and controls, as depicted in Figure 3B, where renal function test and low HDL and hemoglobin displayed a significant impact on differentiating between groups. Regarding the genotype data of four SNPs, RAN (C/T) SNP was the most effective gene variant in ESRD patients. At the same time, controls were mostly affected by DICER1 (A/G) and XPO5 (T/G) variants with a little contributing effect of the DROSHA (A/G) variant (Figure 3C). With the intersection of genes with clinical covariates and laboratory results, clear discrimination between patients and controls was revealed (Figure 3D).

Figure 3

Principal component analysis (PCA) ordination plots for ESRD patient and control. PCA was carried out to demonstrate the similarities and differences among subject groups based on their (A) clinical and demographic characteristics, (B) laboratory results, (C) microRNA machinery gene variants, and (D) integrative data for all variables. Results were plotted on Axes 1 and 2, at which percentage of its component is shown. Dim: dimension. The length of the line indicates its influence and the direction indicates if the influence is positive or negative to the category outcomes. Ellipses are drawn around clusters identified after automatic classification, based on the PCA scores

DISCUSSION

As a part of the miRNA regulome, miRNA processing machinery gene variants might be a source of potential biomarkers. Although several variants in these genes have been explored in other complex diseases, , , , , , virtually no such study has been conducted on their relevance in ESRD.

Here, we identified that the study microRNA processing machinery gene variants, namely, DICER1 rs3742330A/G; RAN rs14035TC/T; and XPO5 rs11077T/G and their combined genotypes might be associated with susceptibility to ESRD under several genetic association models. Furthermore, with the intersection of genes with clinical covariates and lab results, clear discrimination between patients and controls was revealed by multivariate and the principal component analyses. Notably, when analyzing whether the study variants were in HWE, we discovered that most of them were not. Although, it is difficult to speculate the main reason for this, some of the probable causes of population differences shown in the study are selection, relatively small population size, population stratification, and genetic drift. However, the allele frequencies described in our study are comparable to others, as presented in Supplementary Table S2. Since both study groups were not selected from the general population but consisted of hospital‐based individuals chosen as described in the methods section, the above variants were not excluded from analyses. Moreover, the high quality of the genotyping assays (ie, the genotyping call rate =100%, with unambiguous allelic discrimination plots) suggested a violation of HWE assumptions in the study groups rather than technical genotyping errors.

Although the Drosha rs10719*G allele was more prevalent in ESRD patients than the control group, no significant difference was observed in genotype frequencies between the two groups. Interestingly, the presence of this variant was reported to interfere with the interaction between miR‐27b and its target “DROSHA” mRNA. As this variant is located in the miR‐27b‐binding site within 3′UTR of DROSHA, substitution of A by G disrupts the interaction between miR‐27b and DROSHA, increasing the DROSHA expression, which has been associated with increased risk of primary hypertension and preeclampsia. As hypertension is a key player in ESRD etiopathology, with a prevalence of 70%‐80%, it is reasonable to find an increased frequency of the minor allele of this variant in our samples, but not reach the level of significant association as not all enrolled patients had hypertension. This speculation will need further validation on large‐scale cohorts.

The cytoplasmic multidomain DICER1 enzyme is one of the critical regulators of miRNA biogenesis that responsible for pre‐miRNA processing into mature miRNA. Our results revealed that carrying one DICER1 rs3742330*G allele conferred protection against developing ESRD. Although up to the authors’ knowledge, there are no similar published reports, and emerging evidence shows that DICER1 rs3742330 variant may alter its biologic functions and play essential roles in the susceptibility/protection for various diseases and cancers. , , , , , ,

According to in silico analysis, the DICER1 rs3742330 G allele can create new binding sites for three miRNAs (ie, hsamiR‐3622a‐5p, hsa‐miR‐4423‐3p, and hsa‐miR‐5582‐5p), which may trigger degradation of DICER1 mRNA with the corresponding miRNAs. This leads to DICER1 mRNA downregulation in individuals carrying G allele (AG + GG genotypes) as compared to individuals without this allele, which can explain, in part, the significant clinical associations of this variant. Specifically, experimental evidence has investigated the implication of Dicer in atherosclerosis. , It has been reported that in case of reduced endothelial Dicer expression, a decrease of endothelial chemokine ligand 1 expression, monocyte adhesion, atherosclerosis, and the macrophage content of the lesions was observed in “apolipoprotein E knockout mice” after exposure to a high‐fat diet, in part, by reducing the expression of miRNA‐103. Otherwise, Dicer deletion in macrophages (rather than in endothelium) induced advanced atherosclerosis by augmenting beta‐oxidation of fatty acids in foam cells, demonstrating a cell type‐specific differential role of DICER1. Additionally, A allele as part of the haplotype (T‐A‐A) for the DICER1 rs1057035/rs13078/rs3742330 variants was reported to be associated with increased risk of gestational hypertension, while the A‐A haplotype of DICER1 rs13078 and rs3742330 were associated with a protective effect on type 2 diabetes mellitus compared with the A‐T haplotype. All these findings can support, in part, the observed association of DICER1 variant with ESRD development in the study population.

In the present study, RAN rs14035*TT individuals were five times more likely to develop ESRD, while being heterozygote (CT) have twice the risk. Consistent with our finding, recently, it has been reported that patients with type 2 diabetes (which has a significant contribution to the etiopathogenesis of ESRD), who carry TT + TC genotypes, had a 1.89‐fold higher risk of developing macrovascular complications than CC genotype carriers. Ko et al, also demonstrated that RAN rs14035 (CC + CT vs TT) was associated with venous thromboembolism risk in Koreans. As part of the RAN rs14035/rs3803012 haplotype, C‐G was associated with pregnancy‐induced hypertension susceptibility. Given the essential regulatory roles for several cellular processes executed by RAN protein (a member of RAS superfamily of GTPases) and its critical role in the transport of molecules between the nucleus and the cytosol through the nuclear pore complex, in a GTP‐dependent manner, it is not surprised that RAN polymorphism was the most effective gene variant among others with other clinic‐laboratory features in clustering ESRD patients from controls by the multivariate analysis (Figure 3D).

Similar to the protective role of DICER1 G allele, XPO5 (rs11077 T/G) homozygote and heterozygote carriers of G variants were also less likely to develop ESRD in our cases. The previous experiment of knocking down its expression was reported to lead to reduced miRNA levels. Interestingly, Borghini et al, also reported a protective role (with a 32% reduced risk) of this variant against coronary artery diseases in the Italian GENOCOR cohort and they suggested this role is likely to be associated with the circulating levels of vascular miRNAs (ie, miRNA‐132 and miRNA‐140‐3p) and, consequently, their epigenetic function in gene regulation.

Wen et al reported patients with XPO5 rs11077 G allele showed lower XPO5 expression level by luciferase test, and this variant was significantly associated with the onset of thyroid cancer in a Chinese population. Additionally, it has been found that the AA genotype of this variant displayed a trend for high XPO5 expression in esophageal squamous carcinoma tissues by immunochemistry analysis, and these high XPO5 expression levels were also associated with high survival rates among cancer patients. However, in a recent meta‐analysis of 7284 cancer cases, the findings did not support an association between this variant and cancer risk. Due to “the complexity of miRNA biogenesis, the functional biological consequences of these variants remain difficult to determine” as concluded by.

The impact of DROSHA rs10719A/G; DICER1 rs3742330A/G; RAN rs14035C/T; and XPO5 rs11077T/G (A‐A‐T‐T and G‐A‐T‐T) combinations on the increased risk of ESRD was observed in the current population. These 3' UTR gene variants may interfere with mRNA stability by affecting regulatory “protein‐mRNA” and “miRNA‐mRNA” interactions. Accumulating evidence in several diseases, including ones that can impact ESRD directly and indirectly, supported the combined role of these variants and their interactions with certain environmental factors might contribute to disease susceptibility and/or protection in several populations. , , ,

Some limitations of our study merit consideration. Firstly, the study design as a hospital‐based case‐control study can limit the causal relationship. Secondly, other potentially functional variants in miRNA processing machinery genes were not explored with ESRD risk. Additionally, the relatively small sample size, which is reflected in wide confidence intervals observed in particular in association analyses related to combined genotypes, warrants further larger‐scale studies to confirm the findings. Finally, functional studies are required to explore the impact of these variants on gene expression in ESRD cases.

In conclusion, this is the first study showing that the 3'‐UTR DROSHA rs10719, DICER1 rs3742330, RAN rs14035, and XPO5 rs11077 genetic variations and their combined genotypes are associated with ESRD development. These findings highlight the importance of genetic testing among others to be included within the risk stratification panel to help in individualized, targeted therapy and follow‐up of ESRD patients after validation of the findings in larger‐scale, multi‐center, and prospective follow‐up studies in different ethnic populations.

Figure S1‐S2

Click here for additional data file.

Table S1‐S3

Click here for additional data file.