Under the approval of the Ethics Committee of the Nanjing University Medical School-Affiliated Nanjing Drum Tower Hospital, female AIS patients who visited our Joint Scoliosis Center between July 2010 and September 2017 were reviewed for the eligibility to be included in this study. The healthy participants were recruited during the physical examinations for college admission. All the control subjects were excluded to have scoliosis through the Adam’s forward bend test by a senior spine surgeon (YQ). The current case–control association study was composed of 1,980 patients and 2,499 controls from Chinese Han population. Stage 1 was comprised of 980 patients and 1,499 controls that had been recruited in our previous GWAS. The replication stage was comprised of 1,000 patients and 1,000 controls. Informed consent was obtained from the guardians of the participants. Baseline characteristics including initial age and curve magnitude were collected from the medical record.

Genotype imputation was performed with MaCH-Admix software (Liu et al., 2013). The linkage disequilibrium (LD) and haplotype information of the 1000 Genomes Project (phased CHB and CHS data; March 2012 release) were used as the reference, which covered about 40 kb around rs11190870 (Elnitski et al., 2007). After imputation, SNPs with a low imputation quality (R ² < 0.30) or with minor allele frequency (MAF) <0.10 were excluded from the output files. PLINK v1.90 was used to calculate the association of the SNPs with AIS with logistic regression model (Purcell et al., 2007). The online tool LocusZoom was used to plot the genomic regions (Pruim et al., 2010). According to the results of the imputation analysis, the most significant SNP rs1322330 in the promoter region of LBX1 was selected for further replication in 1,000 cases and 1,000 controls.

SNP rs1322330 was genotyped using TaqMan SNP Genotyping Assay. Interpretation of genotyping assay was performed by ABI 7900HT Sequence Detection System (Applied Biosystems, Foster City, CA). Thirty percent of the samples were randomly selected to validate the reliability of the genotyping results.

Blood sample was collected for all participants. Genomic DNA was subsequently extracted using a commercial kit (QIAGEN, Tokyo, Japan) according to the standard protocol. Paraspinal muscles were collected from 134 female AIS patients with main thoracic curve during the surgical intervention. Besides, 28 congenital scoliosis (CS) patients with AIS-like thoracic curve who underwent correction surgery in our clinics were included as the control group. Deep paraspinal muscle biopsies of 1.5 × 1.5 × 1.5 cm³ were collected at the bilateral side of the proximal vertebra of the curve for all the patients. The informed consent of sample collection was obtained from the guardians of all patients.

HEK293T cells were grown in 48-well cell culture plates for 24 h. They were then transiently transfected with a promoterless luciferase vector (pGL4.19-basic) (Promega, United States) or with a pGL4.19-basic vector with the LBX1 promoter fragment (−2,060 to −120) harboring the rs1322330 A-allele or G-allele (construct rs1322330_A-allele or construct rs1322330_G-allele). The cells were transfected with 660 ng of pGL4.19 (with or without insert) along with 33 ng Renilla plasmid. Lipofectamine 2000 (Invitrogen, United States) was used for transfection into HEK293 cells, according to the manufacturer’s protocol. Cells were harvested at 48 h after transfection, and luciferase assays were then performed with the Dual-Luciferase Assay Kit (Promega, United States) according to the manufacturer’s instructions. Cell lysates were tested first for firefly luciferase activity and then for Renilla luciferase activity. Firefly luciferase luminescence values were divided by Renilla luciferase luminescence values from the same transfection to control for differences in transfection efficiency.

Electrophoretic mobility shift assays (EMSA) were conducted using a commercial kit (Viagene Biotech, Co., Changzhou, China) following the manufacturer’s protocol. Nuclear extracts were collected from HEK 293T cells as described previously. We prepared probes for the risk allele G and the non-risk allele A of rs1322330 by annealing 21-bp complementary oligonucleotides and labeling with digoxigenin (DIG)-11-ddUTP (Roche, United States). The DNA/protein binding assay was performed with 10 mg of nuclear extracts using the LightShift Chemiluminescent EMSA Kit (Thermo Fisher Scientific, United States) according to the manufacturer’s recommendations. For competition experiments, nuclear extracts were pre-incubated with excess unlabeled probes. All gel electrophoresis procedures were performed at 4°C. The DNA/protein complexes were detected by streptavidin peroxidase, and signal detection was performed in a 5200 Multi Luminescent Image Analyzer (Tianneng, China).

Twenty AIS patients and 20 CS patients were randomly selected for the muscle fiber analysis. Immunohistochemical staining was performed to determine the ratio of type I fiber to type II fiber of bilateral paraspinal muscles for each patient. Freshly collected muscles tissues were frozen in isopentane/liquid nitrogen. Embedded samples were cryo-sectioned at a thickness of 10 μm for storage at −80°C until further processing. Muscle fiber type distributions were performed on transverse cryosections (10 μm). Serial sections of each muscle were reacted for adenosine triphosphatase (ATPase) at a pH of 9.4. After reacting with the ATPase, there is a clear differentiation into two fiber types. The type I fibers are more lightly stained and the type II fibers more heavily stained. Slides were photographed at 20X (Zeiss, Lab A1), and the cross-sectional area (CSA) was measured for the two fiber types. Approximately 100 fibers were measured per patients within two to three different fields.

Paraspinal muscles of five AIS patients and five CS patients were randomly selected for the isolation of myosatellite cells (MSCs). The MSCs were purified as described previously (Laumonier et al., 2017). After purification, the cells were subsequently cultured in a growth medium with DMEM/F-12 complemented with 20% FBS, 1% α-glutamine, and 1% P/S (penicillin–streptomycin) at 37°C in a 5% CO₂ humidified incubator following the manufacturer’s protocol. Differentiation process was initiated after the cells reached 70%–80% confluence by replacing the growth medium with an equal volume of differentiation medium: DMEM/F-12 complemented with 2% HS, 1% α-glutamine, and 1% P/S. Then, media was changed with a fresh differentiation medium every 2 days. Cell morphology was imaged under the Zeiss Observer A1 system. Cell proliferation assays were performed with a Cell Counting Kit-8 (CCK8). The myogenic differentiation ability was assessed by formation of myotubes, which were immunostained using the myotube-specific MF20 (myosin heavy chain) antibody. The fusion index of myoblasts was evaluated after 7 days and defined as the percentage of nuclei in MF20-stained cells containing two or more nuclei.

Cells were fixed with 4% paraformaldehyde in phosphate-buffered saline (PBS) followed by permeabilization with 0.1% Triton X-100 in PBS. Then, 5% goat serum in PBS was used for blocking. Cells were then incubated with primary antibody (MF20, 1:20, DSHB; Desmin, 1:50, Abcam) overnight at 4°C, followed by incubation with secondary antibody (1:500; Jackson ImmunoResearch) at room temperature for 1 h. Nucleus were counterstained with DAPI solution (1:1,000 dilution in PBS) for 10 min on a rocking platform. Images were taken under 100x with the Leica DM5500 system. Images of stained cells were taken using a fluorescence microscope (Observer A1, Zeiss, Germany).

The MSCs of the control group were transfected with LBX1 lentivirus or with empty vector. The shLBX1-mCherry lentivirus targeting LBX1 and shCtrl-mCherry control lentivirus were obtained from GeneCopoeia (Rockville, United States). The transfected cells were cultured in a growth medium containing 1 μg/ml puromycin (Thermo Fisher Scientific) for another 2 days to select the successfully transfected cells. The sequence of shLBX1 lentivirus targeting LBX1 was as follows: GAC​ATC​CTC​AAC​AAG​CCG​TCT. Knock-down of LBX1 expression was then validated by real-time quantitative PCR (RT-qPCR). All the transfected cells were cultured until the formation of myotube. Cell viability was then compared among the lentivirus-transfected group, the empty vector group, and the blank group.

Freshly collected muscle tissues were snap frozen immediately in liquid nitrogen and stored at −80°C. The total RNA of muscle biopsies and cells was extracted with TRIzol (Thermo Fisher Scientific, United States) according to the manufacturer’s protocol. Reverse transcription of 2 μg of the total RNA was performed with the High Capacity cDNA Reverse Transcription Kit (Applied Biosystems, CA), followed by qPCR with Power SYBR Green (Thermo Fisher Scientific, United States) using Quantstudio^TM 12K Flex real-time PCR platform (Life Technologies, United States). Glyceraldehyde-3-phosphate dehydrogenase (GAPDH) was used as the endogenous control gene for the normalization of mRNA expression. The specific primers are as follows: forward 5′-AGG​ACA​TCC​TCA​ACA​AGC​CG-3′, reverse 5′-CAT​ACC​GTC​GCG​GCC​TTC-3′ for the LBX1 gene; forward 5′-CGG​ACG​TGC​CTT​CTG​AGT​C-3′, reverse 5′-AGC​ACC​TGG​TAT​ATC​GGG​TTG-3′ for the MyoD gene; forward 5′-GGG​GAA​AAC​TAC​CTG​CCT​GTC-3′, reverse 5′-AGG​CGC​TCG​ATG​TAC​TGG​AT-3′ for the MyoG gene; forward 5′-GGA​GCG​CCA​TCA​GCT​ATA​TTG-3′, reverse 5′-ATC​CGC​ACC​CTC​AAG​ATT​TTC-3′ for the MYF6 gene; and forward 5′-GAG​TCA​ACG​GAT​TTG​GTC​GT-3′, reverse 5′-TTG​ATT​TTG​GAG​GGA​TCT​CG-3′ for the GAPDH. All amplifications were performed in triplicate. The mean value of threshold cycle (Ct) scores was calculated to determine the relative expression level. The expression of the target gene was calculated using ΔΔCt method.

Total protein was extracted from the muscle samples or from cultured cells using the lysing buffer (Invitrogen, CA, United States). Protein quantification was conducted using a BCA Protein Assay Kit (Thermo Scientific, CA, United States) following the manufacturer’s instructions. The protein was subjected to SDS-PAGE and transferred onto polyvinylidene difluoride (PVDF) membranes. Sequentially, the membranes were incubated with rabbit anti-LBX1 (ab90836, Abcam), rabbit anti-GAPDH (ab9485, Abcam), and goat anti-rabbit secondary antibodies (ab205718, Abcam) and visualized via enhanced chemiluminescence system. Similarly, rabbit anti-MyoD1 (13812, Cell Signaling Technology), mouse anti-MyoG (F5D, Developmental Studies Hybridoma Bank), and rabbit anti-MYF6 (11754-1-AP, ProteinTech) were used to test the expression of target genes, respectively. GAPDH was used as an internal reference to normalize the quantity of the protein.

The Hardy–Weinberg equilibrium (HWE) test was performed for both patients and controls. The differences of genotype and allele distributions between patients and controls were calculated using chi-square test. The one-way ANOVA test was used to compare the mRNA expression of LBX1 among different genotypes of rs1322330. The Student’s t-test was used to compare the mRNA expression of LBX1 between cases and controls. The Pearson correlation analysis was used to determine the relationship between the tissue expression level of LBX1 and MyoD, as well as the relationship between gene expression and the CSA of myotube in AIS patients. All statistical analysis was conducted with SPSS version 19.0 (SPSS Inc., Chicago, United States). A p value of <0.05 was considered statistically significant.

For the genetic association analysis, the mean age was 15.4 ± 3.7 years (range, 11–18 years) for the patients and 19.1 ± 2.9 years (range, 17–22 years) for the healthy controls. The mean Cobb angle was 35.7° ± 15.2°, ranging from 20° to 72°. For the tissue expression analysis, the mean age was 14.4 ± 2.1 years for AIS patients and 14.3 ± 2.3 years for CS patients (p = 0.82), respectively. For CS patients, 20 had hemivertebra and the other eight were diagnosed as malformation of vertebral body. The mean Cobb angle was 55.1° ± 8.4° for AIS and 57.2° ± 7.7° for CS (p = 0.22), respectively. The proximal vertebral region for sample collection ranged from T3 to T5 level of the thoracic region (Supplementary Table S1).

The imputation analysis indicated that rs1322330 was the most significant SNP located in the promoter region of LBX1 (p = 4.25 × 10⁻⁹) (Figure 1A; Supplementary Table S2). As shown in Table 1, the independent replication analysis confirmed that the frequency of allele A was significantly higher in the patients than in the controls (0.676 vs. 0.618, p = 1.42 × 10⁻⁴). A combination of the imputation analysis and the replication analysis showed that rs1322330 was associated with AIS with genome-wide significance (p = 6.08 × 10⁻¹⁴). Allele A can significantly add to the risk of AIS with an odds ratio (OR) of 1.42 (95% confidential interval: 1.03–1.55).

AIS patients with genotype AA were found to have significantly lower expression of the LBX1 as compared with those with genotype GG (Figure 1C, D). The online databases (HaploReg and UCSC Genome Browser) indicated that rs1322330 was located in a potential TF-binding region with altered signal of DNase cluster (Figure 1B; Kent et al., 2002; Ward and Kellis, 2012). To further address the functional role of rs1322330, reporter gene constructs containing the risk allele (allele A) and non-risk allele (allele G) were prepared. The insertion of the LBX1 promoter fragment consistently resulted in an augmentation of the luciferase activity as compared with the blank construct. The LBX1/A-Luc construct had significantly lower luciferase activity than the LBX1/G-Luc construct (Figure 1E).

The shift band corresponding to allele G probe–protein complexes (Figure 1F, lines 4 and 5) was significantly more intense than that corresponding to allele A probe–protein complexes (Figure 1F, lines 2 and 3), thus suggesting that the two alleles had different affinities for certain transcription factor within the nuclear extracts. A 50-fold excess of unlabeled probes remarkably abrogated the formation of DNA–protein complexes (Figure 1G, lines 3–5), thus confirming the specificity of these interactions.

RT-PCR results are summarized in Supplementary Table S3. AIS patients were found to have remarkably lower mRNA expression of LBX1 in the concave side than in the convex side (Figure 2A). For CS patients, there was no significant difference regarding mRNA expression between the concave side and the convex side. The average mRNA expression of LBX1 was remarkably lower in AIS patients than in CS patients (0.00024 ± 0.000089 vs. 0.00030 ± 0.00012, p = 0.03). Both WB and immunofluorescence analysis confirmed that the protein expression of LBX1 was remarkably lower in AIS muscles than in CS muscles (Figure 2B–F; Supplementary Figure S1).

As shown in Figure 2, there was significantly less type I fiber in the concave muscles as compared with the convex muscles in AIS patients (44.4% ± 12.7% vs. 65.4% ± 16.3%, p = 0.004). By contrast, there was no significant difference regarding the proportion of fiber type between the bilateral sides of paraspinal muscles in CS patients (67.4% ± 21.2% vs. 79.5% ± 25.2%, p = 0.26) (Supplementary Table S4). Intergroup comparison showed that there was globally less type I fiber in the bilateral paraspinal muscles of AIS patients than in CS patients (55.2% ± 14.1% vs. 73.5% ± 22.4%, p = 0.04). In addition, the average CSA of muscle fibers was remarkably smaller in AIS than that in CS (1,350.3 ± 256.1 vs. 1,712.3 ± 374.4 μm², p = 0.02) (Supplementary Table S4). The correlation analysis showed that the mRNA expression of LBX1 was remarkably correlated with the CSA of muscle fibers in the paraspinal muscles of AIS patients (r = 0.375, p = 0.04).

The MSCs were confirmed via the immunostaining of Desmin and Pax7 (Figure 3A–E). The CCK8 assay showed that MSCs isolated from the AIS concave muscles had obviously lower viability than those of CS patients (Supplementary Figure S2A; Supplementary Table S5). Moreover, compared with the CS group, a remarkably lower fusion index of the myotube was observed for the AIS group at the fifth day and the eighth day, respectively (Figure 3F–I; Supplementary Table S6).

To further clarify the involvement of LBX1 in skeletal MSC proliferation and differentiation, we inhibited LBX1 in the MSCs of CS patients with lentivirus, which effectively resulted in the knockdown of LBX 1 by more than 80% (Supplementary Figure S3). A significantly inhibited proliferation rate was observed in the ShLBX1 group by the CCK-8 test (Figure 4; Supplementary Figure S2B; Supplementary Table S7). Moreover, a remarkably decreased formation of the myotube was observed after the knockdown of LBX1 (Figure 4; Supplementary Figure S4).

As the member of myogenesis regulatory family, the expression of MyoG, MyoD, and MYF6 was analyzed at both proliferation and differentiation stages of MSCs. The expression of MyoD was significantly decreased in the ShLBX1 group as compared with the control group (Figure 5A). As for the expression of the other two genes, no significant difference was found between the two groups (Figure 5A). As shown in Supplementary Table S3, the mean expression of MyoD was remarkably lower in the paraspinal muscles of AIS than CS (0.00237 ± 0.00092 vs. 0.00294 ± 0.00129, p = 0.03). There was a significant correlation between the mRNA expression of LBX1 and MyoD (r = 0.57, p = 0.0001) (Figure 5B; Supplementary Table S8). In addition, the expression of MyoD in paraspinal muscles of AIS patients was remarkably correlated with the mean CSA of muscle fibers (r = 0.42, p = 0.04) (Figure 5C).