Identification of genetic variants associated with diabetic kidney disease in multiple Korean cohorts via a genome-wide association study mega-analysis

Background The pathogenesis of diabetic kidney disease (DKD) is complex, involving metabolic and hemodynamic factors. Although DKD has been established as a heritable disorder and several genetic studies have been conducted, the identification of unique genetic variants for DKD is limited by its multiplex classification based on the phenotypes of diabetes mellitus (DM) and chronic kidney disease (CKD). Thus, we aimed to identify the genetic variants related to DKD that differentiate it from type 2 DM and CKD. Methods We conducted a large-scale genome-wide association study mega-analysis, combining Korean multi-cohorts using multinomial logistic regression. A total of 33,879 patients were classified into four groups—normal, DM without CKD, CKD without DM, and DKD—and were further analyzed to identify novel single-nucleotide polymorphisms (SNPs) associated with DKD. Additionally, fine-mapping analysis was conducted to investigate whether the variants of interest contribute to a trait. Conditional analyses adjusting for the effect of type 1 DM (T1D)-associated HLA variants were also performed to remove confounding factors of genetic association with T1D. Moreover, analysis of expression quantitative trait loci (eQTL) was performed using the Genotype-Tissue Expression project. Differentially expressed genes (DEGs) were analyzed using the Gene Expression Omnibus database (GSE30529). The significant eQTL DEGs were used to explore the predicted interaction networks using search tools for the retrieval of interacting genes and proteins. Results We identified three novel SNPs [rs3128852 (P = 8.21×10−25), rs117744700 (P = 8.28×10−10), and rs28366355 (P = 2.04×10−8)] associated with DKD. Moreover, the fine-mapping study validated the causal relationship between rs3128852 and DKD. rs3128852 is an eQTL for TRIM27 in whole blood tissues and HLA-A in adipose-subcutaneous tissues. rs28366355 is an eQTL for HLA-group genes present in most tissues. Conclusions We successfully identified SNPs (rs3128852, rs117744700, and rs28366355) associated with DKD and verified the causal association between rs3128852 and DKD. According to the in silico analysis, TRIM27 and HLA-A can define DKD pathophysiology and are associated with immune response and autophagy. However, further research is necessary to understand the mechanism of immunity and autophagy in the pathophysiology of DKD and to prevent and treat DKD. Supplementary Information The online version contains supplementary material available at 10.1186/s12916-022-02723-4.

necessary to understand the mechanism of immunity and autophagy in the pathophysiology of DKD and to prevent and treat DKD.

Background
Diabetic kidney disease (DKD) is the primary etiology of chronic kidney disease (CKD) in patients with diabetes mellitus (DM) [1] and the leading cause of CKD and end-stage renal disease (ESRD) in most developed countries [2]. Moreover, risk factors for the development of DKD in patients with type 2 DM include cardiovascular risk factors such as high urinary albumin creatinine ratio, old age, hyperglycemia, and hypertension [3]. Although the importance of hyperglycemia in the development of DKD has been illustrated in several studies [4][5][6], some patients with type 2 DM experience a relatively rapid deterioration in renal function, whereas others maintain normal renal function even with suboptimal glycemic levels [7]. Furthermore, patients with DKD showed familial clustering [8] and ethnic group differences [9,10]. This susceptibility highlights the need to identify the specific genetic factors that affect the onset and progression of DKD in patients with DM.
Recently, genome-wide association studies (GWASs) have identified more than 33 genes for DKD in type 2 DM, including APOL1, GABRR1, GCKR, and UMOD [11][12][13][14][15][16][17][18][19]. Moreover, most of the genes reportedly associated with DKD need to be confirmed by further replication studies and a detailed analysis of their functional role in DKD using experimental models [20]. Two complex fundamental features define DKD: the decline of estimated glomerular filtration rate (eGFR) and the presence of proteinuria; hence, it would be better to combine these phenotypes during analysis, to incorporate DM and hypertension-related CKD. However, previous research on GWAS phenotypes for DKD or CKD was confined to the assessment of single phenotypes, such as uric acid, eGFR, ESRD, and proteinuria, and failed to focus on defining DKD [11][12][13][14][15][16][17][18][19]21]. Hypertension is both an underlying risk factor and a consequence of DKD due to persistent high blood pressure in the arteries around the kidney [22]. Additionally, up to 75% of patients with DM also experience hypertension, and individuals with only hypertension frequently exhibit signs of insulin resistance [23,24]. Although previous studies were singularly focused on either CKD or DKD, genes such as UMOD were linked to both hypertensive CKD (non-DKD) and DKD [20,25,26]. Hence, we hypothesized that single-nucleotide polymorphism (SNP)-related traits for DKD could be discovered through a GWAS mega-analysis using multinomial logistic regression (MLR).
In addition, although there are few studies on decreased renal function in middle Eastern descent [21], Japanese [17,27] and Han Chinese [18] populations, most studies were conducted on the European and African-American populations [7,[12][13][14][15][16]19]. There is a need for GWAS on DKD in large-scale Korean multicohorts. As reported in previous studies by the Veterans Health Service Medical Center (VHSMC) [28,29], several elderly veterans are diagnosed with DKD due to the extended duration of type 2 DM. Consequently, studies on DKD in Korean multi-cohorts would be helpful. Toward this goal, we conducted a large-scale GWAS mega-analysis of multi-cohorts, combining the VHSMC cohorts and Korean Genome and Epidemiology Study (KoGES) consortium using MLR with four groups: normal, DM without CKD ("only DM"), CKD without DM ("only CKD") and DKD.

Study population
Clinical and genetic data from multi-cohorts of the VHSMC cohort and the KoGES consortium were integrated in this study (sample size, n = 81,039, Fig. 1). In the previously constructed VHSMC cohort, diagnosed with type 2 DM by VHSMC endocrinologists [28,29], those who met the inclusion criteria (n = 916) were enrolled in this study. The KoGES consortium is a nationwide cohort representative of genome research in Korea [30], of which three cohorts related to population-based studies [Korean Association Resource from Ansan and Ansung (KARE, n = 8840) cohort, KoGES Health Examinees (HEXA, n = 61,568) cohort, and KoGES cardiovascular disease association study (CAVAS, n = 9,715) cohort] were enrolled in this study. This study excluded subjects without Korea Biobank Array genotype data or phenotype (DM, eGFR, albuminuria, and hypertension) data, those who had chronic diseases affecting DM (kidney cancer, pancreatic disease, etc.) and renal function (liver cancer, chemotherapy, etc.), and those younger than 65 years for the control group from the VHSMC cohort and KoGES consortium. After applying the exclusion criteria, 30,069 participants were included in the GWAS mega-analysis (Fig. 1). The VHSMC cohort contained patients diagnosed with type 2 DM and hypertension by certified doctors, whereas in the KoGES cohort data, DM was defined by any of the following four categories (Table 1): (1) DM diagnosis checked in the questionnaire, (2) blood glucose levels ≥ 200 mg/dL 2 h after glucose loading, (3) glycosylated hemoglobin (HbA1c) amount ≥ 6.5%, and (4) overnight fasting blood glucose levels ≥ 126 mg/dL. To develop a distinct control group, participants aged 65 years or above who did not have DM or renal failure were included in the analysis.   [31]. According to the Kidney Disease: Improving Global Outcomes (KDIGO) guidelines, eGFR and albuminuria categories were used to assess a renal complication [32]. Albuminuria was classified into three categories based on the following albumin-to-creatinine ratios: < 3 mg/mmol creatinine, normal to mildly increased; 3-29 mg/mmol, moderately increased; and ≥ 30 mg/mmol creatinine, severely increased. Since the dysfunction of the glomerular barrier (represented by proteinuria) and reduced renal function (assessed using the eGFR) may develop independently, the various phenotypes for DKD were defined as follows (

Genotyping and imputation
Genomic DNA was extracted from venous blood samples, and 100 ng DNA was genotyped using the c Affymetrix Axiom 1.1 (Affymetrix, Santa Clara, CA) [33]. The genotypes were identified using a K-medoid clusteringbased algorithm to minimize the batch effect [34]. The PLINK (version 1.9, Boston, MA) and ONETOOL [35] software packages were used for quality control processes. We excluded samples matching any of the following criteria: (1) sex inconsistencies or (2) a call rate of up to 95%. Furthermore, SNPs were filtered out if (1) the call rate was lower than 95% or (2) the Hardy-Weinberg equilibrium (HWE) test showed P < 1×10 −5 . The genotype imputation was conducted using the Northeast Asian Reference Database imputation server (https:// nard. macro gen. com/), and data of 1779 Northeast Asians [36] were used for the reference panel. Pre-phasing and imputation were performed using Eagle v2.4 [37] and Minimac4 [38], respectively. Post imputation, imputed SNPs were removed if the R-squared value was less than 0.8, there were duplicated SNPs, missing genotype rates were more than 0.05, P-values for HWE were less than 1×10 −5 , or minor allele frequencies were less than 0.01. After quality control and imputation, 6,159,267 SNPs were selected for association analyses.

Genome-wide association analysis
Baseline characteristics of the study population have been reported as means ± standard deviation (SD) for continuous variables and numbers and as proportions for categorical variables. Genome-wide association analyses were conducted with an MLR model for the categorical response variable with four levels (normal, only DM, only CKD, and DKD), implemented in SNPTEST v2.5.6 [39]. We evaluated the overall fit of the model by comparing the likelihood of the two models: a full model with genotype risk factors and a reduced model with covariates only. Age, sex, and ten principal component scores were selected as covariates. To verify that there was no confounding due to population stratification in this study, the variance inflation factor (VIF) was calculated, whereby a VIF value close to 1 indicated no genomic inflation. The regional plot for significant genetic variation was generated using the LocusZoom software with linkage disequilibrium (LD) information of East Asians from the 1000 Genomes Project [40]. The bottom panel displays gene symbols and the location within the region, derived from 1000 genomes (ASN hg19/Nov2014). The threshold for statistical significance in this model was P < 5.0×10 −8 , which is conventionally considered to reflect genome-wide significance. To estimate the relative proportion of phenotypic variance explained by all observed common SNPs, genome-wide complex trait analysis (GCTA v1.91.7) was used for heritability calculation [41].

Conditional analyses
HLA has been reported to have an effect on T1D, and variants in the HLA region were adjusted to remove the confounding effects of T1D on DKD. A conditional MLR model was applied, adjusting the effect of T1Dassociated HLA variants. Among 35 previously reported T1D-associated HLA variants, one SNP (rs9275490) in DR-DQ loci and one SNP (rs9271346) in non-DR-DQ loci were included as covariates for the analysis [42]. Furthermore, to determine the variants that affect DKD independently of DM and CKD, we performed additional logistic regression analyses for DKD after adjusting for the effects of DM and CKD.

Fine-mapping analysis
To prioritize whether the variants discovered are candidate causal variants, fine-mapping analysis was performed. Among significant GWAS hits, even after adjusting for effects of T1D-associated HLA variants, DM, and CKD, only SNPs with genome-wide significantly related to DKD were used for fine-mapping analysis. For each target variant, we first selected the set of SNPs, consisting of the most significant SNP and a 100kb window of SNPs around it. The LD matrix between SNPs was computed using PLINK (version 1.9) [43]. Using a Bayesian approach (PAINTOR method) [44], we estimated the posterior probabilities of causative SNPs at a given fine-mapping locus.

Functional annotation analyses
The eQTL analysis was performed using the Genotype-Tissue Expression (GTEx) dataset. To identify significant eQTL genes, it was assumed that approximately 200,000 SNPs were used in the eQTL analysis considering 20,000 genes and LD block [45]. Therefore, using the Bonferroni correction, we set the significance threshold to 0.05/200,000 (=2.5×10 −7 ). We used the LocusFocus tool to generate a colocalization plot, showing the lead SNP responsible for both GWAS and eQTL signals at loci [40]. GTEx version V7 and 1000 Genomes Phase 3 East Asian LD were used for this plot. The associated genes were further investigated for differently expressed genes (DEGs) in the glomeruli of patients with DM while controlling for age from the Gene Expression Omnibus (GEO) dataset (GSE30529). The platform used for the analysis of GSE30529 was the GPL571 [HG-U133A_2] Affymetrix Human Genome U133A 2.0 Array, which included genes from the kidneys of ten subjects with diabetes and twelve control samples of genes from the kidneys of healthy people. Furthermore, the Search Tool for the Retrieval of Interacting Genes/Proteins (STRING) open-access database was used to identify biological functions based on the identified genes [46]. The minimum required interaction score was set as 0.4 (medium confidence). The Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analysis [47] was also conducted using the STRING database.

Clinical characteristics of the study participants
A total of 33,879 subjects were enrolled in this study, and the baseline characteristics of the study population are presented in Table 2

Genome-wide association analysis for DKD
We conducted a GWAS mega-analysis using MLR on 6,159,267 SNPs from 30,069 subjects from the VHSMC cohort and the KoGES consortium, and they were analyzed after the genotype quality control (Fig. 1). First, we investigated genetic variants that significantly differed between the four groups, using the likelihood ratio test (LRT). Three genetic variants passed the genome-wide significant threshold (P <5.00×10 −8 ) and were identified as novel variations for DKD. The most significant SNP was rs3128852 near OR5V1 (LRT P = 8.21×10 −25 ), followed by rs117744700 near HIATL1 (LRT P = 8.28×10 −10 ), and rs28366355 near human leukocyte antigens HLA-DRB1 and HLA-DQA1 (LRT P = 2.04×10 −8 ; (Table 3 and Additional file 1: Table S1).
Notably, rs3128852 and rs117744700 were significant in DKD even after adjusting for the association of the T1D-associated HLA variants, whereas rs28366355 in the HLA region lost significant association (Additional file 1: Table S2). Furthermore, we confirmed that rs3128852 and rs117744700 were associated with DKD independently of the link between DKD and DM or CKD (Additional file 1: Table S3). Figure 2 represents the quantile-quantile (Q-Q) plot, which verified that there was no inflation in the test statistics (VIF=1.025), and Manhattan plots of the results of the MLR GWAS mega-analysis. Regional association plots are shown in Fig. 3. In addition, two more suggestive variants (P <1.00×10 −6 ) were identified: rs1824125 near PGR (LRT P = 6.58×10 −7 ) and rs75292524 near PARD3B and NRP2 (LRT P = 7.82×10 −7 ; Table 3). The Q-Q and Manhattan plots for six contrasts for each pair are shown in Additional file 2: Fig. S1. Since the HLA region has a notably higher level of variability than the rest of the genome, it was confirmed through the multidimensional scaling (MDS) plot that the false-positive association was not caused by population stratification (Additional file 2: Fig. S2). Furthermore, fine-mapping results supported that the top SNP (rs3128852) with the highest posterior probability (almost 1.0) is likely to have a potential causal effect on DKD (Additional file 1: Table S4 and Additional file 2: Fig. S3A). However, the posterior probability of the second (rs117744700) significant variants was almost zero (Additional file 2: Fig.  S3B). These SNPs may not actually be potential causal variants, or the causality may not be estimated owing to the absence of SNPs with high LD relationships around the target SNP (Fig. 3B, C). Since GWAS is designed for identifying variants associated with the phenotype of interest, and not causal variants, the top SNP can be interpreted as a potential causal SNP for DKD and a related SNP for the other two SNPs. The eQTL  Fig. S4.

Functional annotation
The most significant variant (rs3128852) was an eQTL for TRIM27 (tripartite motif-containing 27) in whole blood (P = 2.10×10 −9 ) and adipose-subcutaneous cells (P = 1.60×10 −7 ; Table 4). However, the second significant variant (rs117744700) did not have any eQTL-associated genes in the GTEx database. The third significant variant (rs28366355) was an eQTL for CYP21A1P, HLA-DQA1, HLA-DQA2, HLA-DQB1, HLA-DQB2, HLA-DRB1, HLA-DRB6, and XXbac-BPG154L12.4 in most tissues. The eQTL analysis results for major tissues that may be associated with renal function (whole blood; artery-aorta, artery-tibial; adipose-subcutaneous; heart-left ventricle; renal-cortex, adrenal gland) are presented in Table 4 and Additional file 1: Table S5. Colocalization plots for the top   Table 5 and Additional file 1: Table S6). Except for HLA-DRB6, which cannot be annotated in the STRING database, five genes were used as input genes in the STRING database to identify known and predicted biological functional networks. The HLA genes (HLA-A, HLA-DQA1, HLA-DQB1, and HLA-DQB2) had a strong interaction with each other (interaction score > 0.75), but TRIM27 did not interact with any of these genes (Additional file 2: Fig. S5). Although 10 KEGG pathways were associated with this network, only two genes (HLA-DQA1 and HLA-A) were annotated for DKD, which are associated with the immune response related to DKD pathogenesis (Table 6). Heritability estimates for DKD, CKD, and DM are presented in Table 7.

Transcriptome-wide association analysis (TWAS) for DKD
mRNA expression levels for 7254 protein-coding genes were imputed for TWAS. Additional file 2: Fig.  S6 presents the Q-Q plot, which verified that there was no inflation in the test statistics (VIF=1.047), and volcano plots of the results of the TWAS. None of the genes showed significant differences between the patients with DKD and the control group in the DEG analysis. Summary statistics for TWAS are provided in Additional file 1: Table S7.

Discussion
In this study, we demonstrated that three novel SNPs (rs3128852, rs117744700, and rs28366355) are significantly linked to DKD. In particular, we noted that the potential causal relationship between rs1328852 and DKD was also confirmed through fine-mapping analysis. The functional analysis for rs3128852 suggests that TRIM27 and HLA-A are potential genes for determining DKD pathophysiology. This study has elucidated the pathological mechanism of DKD through genome analysis. Earlier, GWASs for DKD or CKD were limited to the analysis of specific phenotypes, such as uric acid, eGFR, ESRD, and proteinuria [11][12][13][14][15][16][17][18][19], and several key genes were identified-UMOD, MANBA, DAB2, and SHROOM3 [50]. We hypothesized that comparing the genomes of patients with DKD and healthy normal individuals would reveal DM-related SNPs and CKDrelated SNPs; hence, we divided our investigation into four subgroups, which were specialized for DKD phenotypes. In this study, the eQTL for rs3128852 showed substantial TRIM27 and HLA-A expression, and the results of the subsequent functional genome studies supported these results. LocusZoom plots for significant single-nucleotide polymorphisms. A Regional plot of rs3128852. B Regional plot of rs117744700. C Regional plot of rs28366355. Vertical axis indicates the −log10 of the P-values, whereas the horizontal axis indicates the chromosomal position Each dot represents the single-nucleotide polymorphism (SNP) results for GWAS mega-analysis. Approximate linkage disequilibrium of East Asians from the 1000 Genome Project between the most significant SNPs are listed at the top of each plot; the other SNPs are shown by the 2 key in each plot TRIM27 encodes the tripartite motif protein family, which is involved in a variety of biological activities that may be related to autophagy and pyroptosis [51]. Recent studies demonstrated that TRIM27 was involved in the injury of glomerular endothelial cells in lupus nephritis (LN) through the FoxO1 pathway [52,53] and in IgA nephropathy (IgAN) via T cell signaling [54]. Although the pathogenic processes of DKD and immune-related neuropathy such as LN and IgAN are different, the molecular pathways in cells may overlap, which supports our previous findings that suppression of the protein kinase B pathway could attenuate the damage by mediating the expression of TRIM27 [52,55]. Autophagy is strictly regulated to maintain an optimal balance of cellular component synthesis, degradation, utilization, and recycling of cellular components [56]. When kidney cells are exposed to stress, dysregulated autophagy may contribute to the accumulation of cellular damage, resulting in age-related kidney disease [56]. Several experimental studies have shown that autophagy is inhibited by podocytes or proximal tubule epithelial cells [57][58][59], which is consistent with our results. The accumulation of mitochondria plays a key role in the formation of reactive oxygen species, which activates pro-apoptotic signals and may result in hypertrophy of podocytes [60,61], apoptosis of proximal tubular cells, and kidney fibrosis caused by the WNT-inducible signaling protein-1 [62]. Moreover, upregulation of nephrin expression in the glomeruli inhibits the expression of mammalian target rapamycin, which promotes progressive tubular damage [63]. These results support the notion that profound autophagy dysregulation is related to DKD [64].
Our study discovered that HLA-A-related genes (HLA-A, HLA-DQA1, HLA-DQB1, HLA-DQB2, and HLA-DRB6) were involved in the etiology of DKD. Although research has been conducted on the role of the immune system in CKD development, there are few studies on DKD; therefore, these results must be taken into account. According to previous studies [65,66], renal function is associated with HLA type, such as HLA-A*01:01, HLA-A*03:01, and DQB1*02:01, which were related to CKD or ESRD. Furthermore, immune mechanisms may play a crucial role in DKD pathogenesis, especially leukocyte accumulation and associated molecular mechanisms [67]. Moreover, hyperglycemia-induced oxidative stress pathologically stimulates circulating immune cells, which enter the affected kidney and exacerbate tissue inflammation by producing pro-inflammatory cytokines and chemokines abundantly [68]. Furthermore, DNA methylation via the upregulated activity of DNA methyltransferase 1 revealed that inflamed memory immune cells aggravate DKD [54]. In addition, HLA-DPA1 may be involved in immune mechanisms underlying DKD development, since it has been identified as a significant gene for DKD [69]. Nevertheless, inflammatory response is a major factor in the progression of DKD, and the immune response exacerbates inflammation, indicating that the adaptive immune response is crucial in DKD [70,71]. Combining the relevance of immune responses in DKD and the results of this study, immune responses and autophagy may be considered as possible pathways in the pathophysiology of DKD.
The true strength of this study includes the utilization of a relatively large elderly cohort sample that provides a better DKD phenotype. This is because diabetes is an age-related disorder, and a longer duration of diabetes   is connected with an increased risk of developing DKD [24]. Furthermore, our study focused on DKD using MLR and compared it with DM without CKD and CKD without DM, which have not been evaluated earlier. However, there are a few limitations of this study. First, in this study, DKD was not diagnosed through kidney biopsy, but by clinical diagnosis. In a previous genetic study of patients with DKD selected based on a definite diagnosis via a kidney biopsy [72], it was difficult to recruit a significant number of participants due to the diversity of the etiological processes of DKD, limiting the study methodology. Moreover, diagnostic kidney biopsy is rarely performed in clinical practice [73]. However, our work has overcome this limitation by using in silico analysis and generating reproducible results, minimizing this disadvantage. Second, subjects with T1D would be at greater risk of developing DKD due to longer DM duration, and the genome-wide significant SNPs can have a vertical or horizontal pleiotropic effect on T1D and DKD. For instance, rs28366355, located near HLA, was significantly associated with DKD in our analysis; this significant result can be inferred from its association with T1D. However, for rs3128852 and rs117744700 or nearby SNPs, no significant results were reported and pleiotropic effects thereof on T1D and DKD are not expected. Furthermore, the prevalence of T1D is very rare, ranging from 0.017 to 0.021% in Koreans [74]. Thus, most individuals with DKD in our analyses may not experience T1D, and there is very low chance of a confounding effect by T1D. Further studies with individuals with DKD and T1D and T2D disease status are necessary. Third, as a study design, a mega-analysis was conducted; however, the disparity in the cohort mix is limiting. Since the CAVAS cohort was established for cardiovascular disease research, it has relatively high prevalence of hypertension and non-DM-CKD. Furthermore, the VHSMC cohort is a hospital-based cohort, whereas the KoGES consortium is based on community survey results. Hence, the severity of diabetes in these cohorts is different. By contrast, the difference in diabetes severity across these cohorts might yield more relevant results when analyzed with respect to real-world data. Fourth, bioinformatics analysis revealed that certain genetic variants and metabolic pathways were related to DKD pathogenesis, but the underlying mechanism of these factors needs further investigation. Because a simple overlap of GWAS lead variants with GTEx nominal P-value results is expected to yield several false-positive candidate causal genes [75], we conducted a QTL colocalization analysis and TWAS. However, the results were underpowered given the sample size in our study. In the PheWAS catalog (https:// phewa scata log. org/), for the SNP located in HLA-A (rs2860580) and its LD relationship with the top SNP (rs3128852), there is a significant association with genitourinary phenotype (Additional file 2: Fig. S7). Further studies are required to elucidate the mechanism through which immune responses and autophagy influence DKD pathogenesis, as discovered in our study. Fifth, our study results are related to some HLA-related regions, and typically, the HLA region has a higher level of variability than the rest of the genome; thus, we need to be careful when interpreting the results. To address these concerns, we attempted to address this issue by applying an MDS plot and conducting a fine-mapping analysis. Nevertheless, the variants identified as significant genome-wide in our study were not all fine-mapped because fine-mapping analysis requires high-quality genetic data and a much larger sample size than that required for a GWAS [76].

Conclusions
This study has demonstrated that three novel SNPs (rs3128852, rs117744700, and rs28366355) are significantly associated with DKD based on the MLR GWAS mega-analysis. Moreover, the causal relationship between rs1328852 and DKD was confirmed through fine-mapping analysis. The functional analysis of the genetic variants (rs1328852) detected has revealed that TRIM27 and HLA-A, associated with immune response and autophagy, contribute to the etiology of DKD. Considering the mechanism through which immune responses and autophagy influence the pathophysiology of DKD, further research is necessary for effective prevention and treatment of DKD.
Additional file 1: Table S1. Results of the GWAS mega-analysis using multinomial logistic regression (P < 0.00001). Table S2. Results from multinomial logistic regression with participant category as the outcome and adjusting for known T1D-associated HLA markers (rs9275490 and rs9271346) [42]. Table S3. Results of the logistic regression for the top two SNPs (rs3128852 and rs117744700) after adjusting for the effects of DM and CKD. Table S4. Results of fine-mapping analysis using top 3 SNPs .  Table S5. Results of eQTL analysis using top 3 SNPs from GTEx. Table S6.
Results of DEG analysis. Table S7. TWAS results for DKD.

Acknowledgements
This study was conducted using the bioresources from the National Biobank of Korea, the Korea Disease Control and Prevention Agency, Republic of Korea (KBN-2021-042). Korean Genome and Epidemiology Study (KoGES) consortium was used to obtain the following three population-based cohorts.