Genome-Wide Association Study Adjusted for Occupational and Environmental Factors for Bladder Cancer Susceptibility

This study examined the effects of single-nucleotide polymorphisms (SNPs) on the development of bladder cancer, adding longest-held occupational and industrial history as regulators. The genome purified from blood was genotyped, followed by SNP imputation. In the genome-wide association study (GWAS), several patterns of industrial/occupational classifications were added to logistic regression models. The association test between bladder cancer development and the calculated genetic score for each gene region was evaluated (gene-wise analysis). In the GWAS and gene-wise analysis, the gliomedin gene satisfied both suggestive association levels of 10−5 in the GWAS and 10−4 in the gene-wise analysis for male bladder cancer. The expression of the gliomedin protein in the nucleus of bladder cancer cells decreased in cancers with a tendency to infiltrate and those with strong cell atypia. It is hypothesized that gliomedin is involved in the development of bladder cancer.


Introduction
There were approximately 83,730 new cases of bladder cancer (64,280 in men and 19,450 in women) and approximately 17,200 deaths from bladder cancer (12,260 in men and 4940 in women) in the United States in 2021. The rates of new bladder cancer and death due to bladder cancer have been decreasing slightly in women in recent years, whereas in men, incidence rates have been decreasing, but death rates are stable [1]. In Japan, 23,230 cases (17,555 in men and 15,675 in women) of bladder cancer were diagnosed in 2018, and the number of deaths from bladder cancer was 8911 (6014 in men and 2897 in women) in 2019 [2].
Muscle invasive bladder cancer requires highly invasive treatments such as radical cystectomy and systemic chemotherapy. In addition, even non-muscle infiltrating bladder cancer often recurs in the bladder and requires multiple treatments. Thus, medical treatment for bladder cancer requires a great deal of time and medical expenses.
It is well known that smoking is a risk factor for developing bladder cancer [3,4]. Regarding alcohol drinking, the American Society of Clinical Oncology stated in 2018 that more than 5% of new cancer cases were due to alcohol consumption [5]. We also reported that alcohol consumption is an independent risk factor for the development of bladder cancer in the Japanese population [6].
Occupational and environmental factors are important, in addition to genetic predisposition for bladder cancer. There was a 45-year observational study in Nordic countries on the association between occupation and the development of bladder cancer. According Genes 2022, 13, 448 2 of 16 to this study, occupations with a significantly increased incidence of urothelial cancer, with a standardized incidence ratio of 1.20 or higher, include male waiters, chimney sweeps, hairdressers, assistant nurses, seamen, plumbers, cooks and stewards, beverage workers, female tobacco workers, printers, waiters, chemical process workers, sales agents, hairdressers, mechanics, and administrators [7].
This study examined the effects of single-nucleotide polymorphisms (SNPs) in the germline genome on the development of bladder cancer in Japan, adding occupational and industrial history as regulators.
A genome-wide association study (GWAS) comprehensively searches the entire genome for gene polymorphisms that exhibit significant frequency differences between an unrelated patient population of a specific disease and an unrelated control population.
In genome-wide studies that analyzed the genomes of bladder cancer patients, 57 SNPs that may increase susceptibility to bladder cancer were identified in the GWAS Catalog (Supplementary Table S1). In addition, there are many GWAS papers on bladder cancer [8][9][10][11][12][13][14][15][16][17][18][19][20][21][22][23][24][25]. In particular, the NAT2 slow acetylator and GSTM1 null genotype are considered to be potential genetic risk factors for the development of bladder cancer [8]. Polymorphisms in the NAT2 gene were also investigated in Japan, with a risk ratio of 7.80-times [26]. In addition, a relatively large GWAS for Japanese bladder cancer patients was announced in 2015, and although smoking has been examined and adjusted as an environmental factor, occupational factors have not been examined [21]. Therefore, it is important to examine the relationship between bladder cancer in the Japanese population and SNPs by adjusting for the industrial/occupational history, in addition to sex, smoking history, and alcohol drinking history.

Materials and Methods
The genome was purified using 10 mL of blood mixed with EDTA collected from 352 bladder cancer patients (302 males, 50 females) and 434 control patients (395 males, 39 females) at Japan Organization of Occupational Health and Safety, Kanto Rosai Hospital and Tokyo Metropolitan Tama Medical Center. Control patients did not include those with upper tract urothelial cancer because bladder cancer and upper tract urothelial cancer are considered to be malignant tumors that are anatomically, histologically, and epidemiologically similar.
Occupational and environmental data were obtained from the Inpatient Clinico-Occupational Database of Rosai Hospital Group (ICOD-R), provided by the Japan Organization of Occupational Health and Safety. The ICOD-R includes an occupational history of current and past three jobs, information on smoking, and alcohol habits using interviews and questionnaires completed at the time of admission. Detailed occupational histories were coded with three-digit codes in the Japan Standard Occupational Classification and Japan Standard Industrial Classification corresponding to the International Standard Industrial Classification and International Standard Occupational Classification, respectively [27]. The Japan Standard Occupational Classification is composed of 12 major groups, 74 minor groups, and 329 unit groups [28], whereas the Japan Standard Industrial Classification is composed of 20 divisions, 99 major groups, 530 groups, and 1460 industries [29]. Other clinical data were obtained from electronic medical records. Missing values exist due to omission or lack of description by patients.

A New Classification of Industry/Occupation
To create a new classification, we divided the occupations into four groups: professional, service, management, and blue-collar workers, and further divided the industries into three groups: white-collar industry, blue-collar industry, and service industry. These two kinds of groups were combined into a total of 12 (4 × 3) industry/occupation classes [30,31]. Using this classification, tentatively named the Zaitsu classification, we previously reported an association between occupation and the prognosis of bladder cancer [32].

Clinical and Environmental Factor
From the clinical data, categorical variables were preliminarily analyzed by Fisher's exact test between two or multiple groups, and continuous variables were preliminarily analyzed by the Mann-Whitney U test. Furthermore, logistic regression analysis was performed with the development of bladder cancer as the objective variable, whereas age, sex, Brinkman index (BI) classified into four groups (0: BI 0, 1: 1-399, 2: 400-799, 3: 800≤), alcohol consumption history (2 levels, yes or no), and industrial / occupational classifications of the longest-held job for each patient were explanatory variables.
The  Table S3). From the logistic regression models of a, b, and c, the explanatory variables related to industry/occupation were selected by the backward step-wise method using the Akaike information criterion.

Genotyping
Performed by Riken Genesis Co., Ltd. (Taito-Ku, Tokyo, Japan). Samples were genotyped using the Illumina Infinium Asian Screening Array-24 v1.0 BeadChip, which combines genome-wide coverage of East Asian populations, relevant clinical research content, and scalability for genomic screening. For quality control of samples, we excluded those with (i) a sample call rate < 0.99, (ii) a person with the lowest call rate from the pairs with a proportion IBD (identity-by-descent) > 0.1875, and (iii) outliers from Japanese clusters identified by principal component analysis using the genotyped samples and East Asians in the International Genome Sample Resource [33] (The 1000 Genomes Project Consortium 2015). For quality control of genotypes, we excluded those with a (i) SNP call rate < 0.99 or (ii) p-value for the Hardy-Weinberg Disequilibrium test < 0.001.

Imputation
We utilized SNP imputation for all samples under 1000 Genomes Project Phase 3 as a reference panel [34]. We implemented the pre-phasing by Eagle [35,36] and imputation by Minimac3 [37]. After imputation, we excluded SNPs with an imputation quality of R-square < 0.3.

GWAS
We conducted 6 GWAS patterns for bladder cancer development by logistic linear models using SNP dosage obtained by SNP imputation and Efficient and Parallelizable Association Container Toolbox (EPACTS) [38]. In the association test, age, sex, smoking history (Brinkman Index, ordered category with 0 < 1 < 2 < 3 levels), alcohol consumption history (2 levels, yes or no), and several patterns of industrial/occupational classifications were added to logistic regression models. Tested industrial/occupational classifications were: (i) 1 variable with 20 levels for industrial classification divisions, (ii) selected industrial classification division(s) from 20 variables with 2 levels (yes or no) by the backward step-wise method in a logistic regression model without SNP dosage, (iii) 1 variable with 12 levels for occupational classification major groups, (iv) selected occupational classification major group(s) from 12 variables with 2 levels (yes or no) by the backward step-wise method in a logistic regression model without SNP dosage, (v) selected industrial classification major groups from 35 variables with 2 levels (yes or no) by the backward step-wise method in a logistic regression model without SNP dosage, and (vi) the Zaitsu classification. We also used only male samples for GWAS, taking into account sex differences in some occupations. We did not conduct GWAS using only females due to the small number of cases. We set the genome-wide significance level for our study at p = 5 × 10 −8 and suggestive association level at p = 10 −5 [39].

Gene-Wise Analysis
For SNPs contained within 50 bp upstream and downstream of the gene regions defined in Ref Gene [40], we calculated the genetic score (GS) [41] as described below, and the association test between bladder cancer development and GS for each gene region was evaluated by the Burden test and SKAT-O test using EPACTS (version 3.2.6) (University of Michigan, Ann Arbor, MI, USA). In our study, we performed gene-wise analysis for 20,865 regions. We set the genome-wide significance level for our study at p = 2.4 × 10 −6 (=0.05/20,865) and suggestive association level at p = 10 −4 . Adjusting factors in GWAS were also included in the gene-wise analysis. GWAS and gene-wise analysis were performed by StaGen Co., Ltd. (Taito-ku, Tokyo, Japan 111-0051).
Here, the GS i of an individual patient is equal to the weighted sum of the individual's genotypes, x j (0, 1, 2), at SNPs in gene i . Weights (β j ) are calculated by EPACTS and M is the number of SNPs in gene i .

Immunohistochemistry
The expression of gliomedin protein was examined by tissue immunostaining using paraffin-embedded bladder tumor tissue removed by transurethral resection of the bladder tumor. The antibody used was anti-GLDN (gliomedin) polyclonal antibody (26185-1-AP, Proteintech, Rosemont, IL, USA). Two independent pathologists evaluated histological staining by the immunoreactive score [42] and individual scores were analyzed after averaging.

Study Approval
The Ethical Committee of the Japan Organization of Occupational Health and Safety approved the experiments (2018-2). All experiments were performed in accordance with relevant guidelines and regulations, including any relevant details. Written informed consent was received from patients prior to inclusion in the study.

Clinical and Environmental Factors
The age of the bladder cancer patients included was slightly lower than the control patients for men and higher for women (Table 1). Malignant tumors other than urothelial cancer were observed in 13.9% of men and 18.0% of women in the bladder cancer group, and 72.2% of men and 59.0% of women in the control group. In the male control group, 46.6% had prostate cancer and 11.1% had renal cell carcinoma (Supplementary Table S4).
In terms of smoking history, a high Brinkman index classified into four stages and the development of male bladder cancer were slightly related. The Brinkman index 2-3 group had more male bladder cancer than the Brinkman index 0-1 group ( Table 1). As for alcohol consumption history, bladder cancer patients drank slightly less alcohol than control patients overall ( Table 1). The overall distributions of the divisions of industrial classification (Table 2), occupational classification major groups (Table 3), and groups in the Zaitsu classification (Table 4) were not significantly different from controls in male, female, and all bladder cancer patients. Looking at the individual divisions of industrial classification, bladder cancer was less frequent in division G and more frequent in division S in male cases and all cases ( Table 2). In addition, in the individual major groups of occupational classification, bladder cancer was significantly more common in the major group F in male cases and all cases (Table 3).        In the selection of explanatory variables concerning industry/occupation by the backward step-wise method for logistic regression models, industrial classification divisions G and S, and occupational classification major group F, were selected in male cases from all divisions and from all major groups, respectively, whereas in all cases including both men and women, industrial classification divisions G, L, and S, and occupational classification major group F remained.
In addition, from the major groups in industrial classification divisions D, E, and H, "Manufacture of general-purpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport" were selected as explanatory variables in male cases, whereas "Miscellaneous manufacturing industries", "Construction work by specialist contractor", and "Railway transport" remained in all cases. The overall distribution of bladder cancer cases in these industrial major groups was not different from that of the controls (Table 5). Table 5. Distribution of industrial classification major groups in divisions D, E, and H in males.

Sample QC
For 789 genotyped samples of 830 samples, 21 were excluded in which the sample call rate was <0.99, the proportion IBD was >0.1875, and outliers from Japanese clusters identified by principal component analysis (Supplemental Figure S1). In our study, we used 766 samples in GWAS (Supplementary Table S5).

SNP QC for Imputation
We selected SNPs to be used for SNP imputation. The number of SNPs loaded on the chip was 659,184 and the number of SNPs genotyped was 657,060. In addition, there were 641,043 SNPs with a definite chromosomal location, 395,708 SNPs with a call rate of 99% or higher, a p-value of 0.0001 for the Hardy-Weinberg law of equilibrium, and a minor allele frequency (MAF) of 1% or higher. SNP imputation was performed using 395,708 SNPs (Supplementary Table S6).

Imputation
The number of SNPs able to be imputed was 47,109,297. Of the 47,109,297 SNPs, 11,175,945 with an R-squared value greater than 0.3 were used in the GWAS (Supplementary Table S7).

Results of GWAS and Gene-Wise Analysis
No SNPs satisfying the genome-wide significance level 5 × 10 −8 were detected in GWAS in all cases or in male cases. A Manhattan plot of the genome-wide association test for each analysis pattern is shown in Figure 1. In the gene-wise analysis, no genes satisfying the genome-wide significance level 2.4 × 10 −6 were detected in all cases or in male cases. In GWAS and gene-wise analysis, the gliomedin (GLDN) gene located at 15q21.2 satisfied both a suggestive association level of 10 −5 in GWAS and suggestive association level of 10 −4 in gene-wise analysis (Table 6 and Figure 2).

Imputation
The number of SNPs able to be imputed was 47,109,297. Of the 47,109,297 SNPs, 11,175,945 with an R-squared value greater than 0.3 were used in the GWAS (Supplementary Table S7).

Results of GWAS and Gene-wise Analysis
No SNPs satisfying the genome-wide significance level 5 × 10 −8 were detected in GWAS in all cases or in male cases. A Manhattan plot of the genome-wide association test for each analysis pattern is shown in Figure 1. In the gene-wise analysis, no genes satisfying the genome-wide significance level 2.4 × 10 −6 were detected in all cases or in male cases. In GWAS and gene-wise analysis, the gliomedin (GLDN) gene located at 15q21.2 satisfied both a suggestive association level of 10 −5 in GWAS and suggestive association level of 10 −4 in gene-wise analysis (Table 6 and Figure 2). (ii) selected industrial classification divisions G and S for male bladder cancer, with divisions G, L, and S for all bladder cancer; (iii) 1 variable with 12 levels for occupational classification major groups; (iv) selected occupational classification major group F; (v) selected industrial classification major groups in D, E, and H, i.e., "Manufacture of general-purpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport" for male bladder cancer, and "Miscellaneous manufacturing industries", "Construction work by specialist contractor", and "Railway transport" for all bladder cancer; and (vi) the Zaitsu classification, The SNPs in the GDLN region are plotted in red.  (ii) selected industrial classification divisions G and S for male bladder cancer, with divisions G, L, and S for all bladder cancer; (iii) 1 variable with 12 levels for occupational classification major groups; (iv) selected occupational classification major group F; (v) selected industrial classification major groups in D, E, and H, i.e., "Manufacture of general-purpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport" for male bladder cancer, and "Miscellaneous manufacturing industries", "Construction work by specialist contractor", and "Railway transport" for all bladder cancer; and (vi) the Zaitsu classification, The SNPs in the GDLN region are plotted in red. Table 6. Results of GWAS and gene-wise analysis for bladder cancer. Results that satisfied both p < 10 −5 by GWAS and p < 10 −4 by gene-wise analysis were selected. Analysis: industrial/occupational factors added in GWAS and gene-wise analysis in male bladder cancer: (i) 1 variable with 20 levels for industrial classification divisions; (ii) selected industrial classification divisions G and S for male bladder cancer; (iii) 1 variable with 12 levels for occupational classification major groups; (iv) selected occupational classification major group F; (v) selected industrial classification major groups in D, E, and H, i.e., "Manufacture of general-purpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport", Chr: chromosome, BP: base pair position, Ref: reference allele, Alt: alternative allele, SE: standard error. tional classification major groups; (iv) selected occupational classification major group F; (v) selected industrial classification major groups in D, E, and H, i.e., "Manufacture of general-purpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport", Chr: chromosome, BP: base pair position, Ref: reference allele, Alt: alternative allele, SE: standard error. Figure 2. A representative regional plot of GLDN region. Added industrial/occupational factors were selected among industrial classification major groups in D, E, and H, i.e., "Manufacture of general-purpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport" for male bladder cancer.
In the Manhattan plot of male bladder cancer cases, there were peaks satisfying a suggestive level of p < 10 −5 between LINC00922 and CDH5 in chromosome 16, and between LINC00473 and PDE10A in chromosome 6, but they were outside the genetic regions in the regional plots of GWAS results (Supplementary Figure S2).
In all bladder cancer cases, Manhattan plot peaks satisfying p < 10 −5 were observed near WNT2B in chromosome 1 and XYLB of chromosome 3 (Supplementary Figures S3  and S4), but they did not satisfy the suggestive association level of p < 10 −4 in gene-wise analysis.

Expression of Gliomedin Protein in Bladder Cancer Tissues
The expression of the gliomedin protein (Figure 3) in the nucleus of bladder cancer cells was lower in cancers with a tendency to infiltrate and those with strong cell atypia, as shown in Table 7. The expression of the gliomedin protein in the cytoplasm of cancer cells and in the nucleus of stromal cells was not associated with the degree of cancer infiltration or cancer cell atypia. Figure 2. A representative regional plot of GLDN region. Added industrial/occupational factors were selected among industrial classification major groups in D, E, and H, i.e., "Manufacture of generalpurpose machinery", "Miscellaneous manufacturing industries", "Construction work by specialist contractor", "Equipment installation work", and "Railway transport" for male bladder cancer.
In the Manhattan plot of male bladder cancer cases, there were peaks satisfying a suggestive level of p < 10 −5 between LINC00922 and CDH5 in chromosome 16, and between LINC00473 and PDE10A in chromosome 6, but they were outside the genetic regions in the regional plots of GWAS results (Supplementary Figure S2).
In all bladder cancer cases, Manhattan plot peaks satisfying p < 10 −5 were observed near WNT2B in chromosome 1 and XYLB of chromosome 3 (Supplementary Figures S3 and S4), but they did not satisfy the suggestive association level of p < 10 −4 in gene-wise analysis.

Expression of Gliomedin Protein in Bladder Cancer Tissues
The expression of the gliomedin protein (Figure 3) in the nucleus of bladder cancer cells was lower in cancers with a tendency to infiltrate and those with strong cell atypia, as shown in Table 7. The expression of the gliomedin protein in the cytoplasm of cancer cells and in the nucleus of stromal cells was not associated with the degree of cancer infiltration or cancer cell atypia. Genes 2022, 13, x FOR PEER REVIEW 11 of 16

Discussion
Kawasaki City, where Kanto Rosai Hospital is located, is a traditional heavy industry area adjacent to the Tokyo Metropolitan area. Fuchu City, where Tokyo Metropolitan Tama Medical Center is located, is a commercial and residential area on the outskirts of Tokyo. Therefore, workers engaged in the primary sector of industry and mining industry were limited in this study.
The reason why the industrial classification divisions G, L, and S, and the occupational classification major group F were particularly adopted as the adjusting factors in certain patterns of GWAS analysis is that these factors were selected by the backward step-wise method in the analysis of the patient background. In addition, the adjusting factors were selected from the industrial classification major groups included in the industrial classification divisions D, E, and H for men because a relatively large number of cases were included in these three divisions.
The occupations vulnerable to bladder cancer in previous reports were fairly specific and limited. Compared with these, this study mainly used the relatively rough classification of industry/occupation such as industrial classification divisions, major groups, and occupational classification major groups. Therefore, the purpose of this study was not to examine the relationship between SNPs and specific environmentally exposed substances, such as nicotine and aromatic amines, but rather to incorporate the contribution of broader industrial/occupational environmental factors, such as stress stimulation and work environment, into the development of bladder cancer as adjustment factors. Under these conditions, the gliomedin gene was detected in this study by GWAS and gene-wise analysis as a gene that may be associated with the development of bladder cancer in males.
GWAS is widely performed to replicate obtained results with other datasets. However, detailed recording of occupational/industrial history, such as ICOD-R, is not comprehensively enforced in Japan, making it difficult to replicate GWAS with occupational/industrial history as an adjusting factor. Therefore, in this study, in addition to GWAS to verify one SNP, the results were supported by performing gene-wise analysis to examine the association between a given pathological condition and a certain gene as a whole.
Kaneko et al. used ICOD-R occupational classification major groups to demonstrate that occupations with high physical activity reduced the risk of cancer [27]. They also compared the categories included in the manufacturing industry division (Division E) of ICOD-R and noted that the incidence of ureter cancer in the electronics category is higher than that in the food manufacturing category [43]. Therefore, adding the industrial/occupational classification to the adjusting factors of GWAS, even if it is relatively rough, is considered to be meaningful in examining the development of cancer.
The control group in this analysis included several malignant tumor diseases other than urothelial cancer. The inclusion of many cases of other malignancies in the control group of GWAS for bladder cancer is controversial. It is thought that pathways common to malignancies in general are less likely to be detected, but on the other hand, it may be more effective in order for pathways specific to bladder cancer to emerge.
The gliomedin gene encodes a protein containing olfactomedin-like and collagen-like domains. The gliomedin protein, which is present in both transmembrane and secretory forms, promotes the formation of the Node of Ranvier in the peripheral nervous system [44]. Mutations in the gliomedin gene cause lethal congenital contracture syndrome [45]. Autoantibodies to the gliomedin protein have also been identified in patients with multifocal motor neuropathy serotypes [46]. An important paralog of the gliomedin gene is olfactomedin protein family [47].
The expression of gliomedin mRNA and protein is found in the nuclei of many types of cancer cells, including urothelial cancer [48]. The early deregulation of gliomedin during liver tumorigenesis was previously reported [49]. The gliomedin paralog, olfactomedin 4, is a glycoprotein with an olfactomedin-domain, which is involved in numerous intracellular signaling pathways, including NF-κB, and is associated with innate immunity. Furthermore, olfactomedin 4 suppresses the development and progression of cancer [50,51], and tumorigenesis is observed in olfactomedin4 deficient mice [52]. Thus, it is hypothesized that gliomedin is also involved in the development of bladder cancer by a mechanism similar to that of olfactomedin 4 in innate immunity and oncogenesis in a certain environment.
As the expression of the gliomedin protein in the nucleus of cancer cells is decreased in bladder cancer with strong nuclear atypia and infiltration tendency, it is speculated that gliomedin may act as a tumor suppressor factor in bladder cancer. As GWAS suggested a relationship between the gliomedin gene and the development of bladder cancer in men, further studies are required.

Conclusions
In conclusion, the gliomedin gene was suggested to be related to the development of male bladder cancer by adding longest-held occupational and industrial history as regulators in the GWAS and gene-wise analysis. In addition, the expression of the gliomedin protein in the nucleus of bladder cancer cells was lower in cancers with a tendency to infiltrate and those with strong cell atypia.