The proband was recruited and examined at the Department of Neurology of the First Affiliated Hospital of Fujian Medical University. Detailed neurological examinations were performed on all members of the family by at least two senior neurologists. Blood samples were collected according to standard methods. Written informed consent was obtained from all participants in this study. The research project was approved by the Ethics Committee of the First Affiliated Hospital of Fujian Medical University.

Genomic DNA was extracted from the peripheral blood cells of patients and normal controls using Qiagen kits. Whole exome DNA was captured using a SureSelect Human All Exon V6 kit (Agilent) and sequenced using the Illumina HiSeq 3000 platform as described in our previous study (17). Fragment sequences were aligned to the consensus sequence (UCSC hg38). Variant calling was performed by Genome Analysis Toolkit (GATK) (18) and annotated with ANNOVAR (19). We filtered the SNP sets obtained from the total sequencing data using the parameters of GQ>20 and DP>8 for genotype quality (GQ) and read depth (DP), respectively. Variants were excluded based on the following criteria: (a) the variants did not affect the amino acid sequence; (b) the allele frequency was >1% in the 1,000 Genomes Project, ESP database, or gnomAD. Then we screened potential disease-causing genes according to their patterns of inheritance.

Total RNA was isolated from the peripheral blood leukocytes of family members. RNA was extracted using a Trizol extraction kit (Invitrogen, Carlsbad, CA, USA) and then synthesized to cDNA with a PrimeScript RT reagent kit (TAKARA BIO, Kusatsu, Shiga, Japan) according to the manufacturer's protocols. Primers were designed to target ALDH18A1 exons 1–3 to confirm aberrant splicing (forward:5′-ACGGAAGAAAAAAGAGAGTGAG-3′; reverse: 5′-AAGGACTTGCCATGTGTACG-3′). The products amplified by the above primers were separated by agarose gel electrophoresis. The DNA fragments were then purified from the gels and sent for Sanger sequencing.

The products amplified by the above primers were used for high-throughput sequencing (HTS). The detailed methods were described in the previous literature (20). In brief, the sequencing library was constructed (New England Biolabs), quality checked, and sequenced with paired-end 150-bp reads using the Illumina HiSeq X-Ten platform (Novogene). Paired end reads of 150 bp were pre-processed to remove adaptors and low quality reads using fastp (v 0.19.6) (21). The following criteria were used to remove the low-quality reads: (i) reads containing more than 15 N bases, (ii) a quality score of < 20, and (iii) more than 50% of the sequenced bases have low quality scores. Clean reads were subsequently aligned to the appropriate reference sequence using the Smith-Waterman algorithm. Mutation types were extracted from the aligned results and the number of all types of mutations were also counted. Subsequently, percentage of diverse transcripts was calculated as the number of reads align to mutant sequence over the total number of reads covering the same locus.

The blood samples from the patient carrying ALDH18A1 mutations and three gender matched healthy controls were centrifuged and the plasma was fractionated prior to storage at −80°C. The concentration of P5CS, encoded by ALDH18A1, was measured by human enzyme linked immunosorbent assay (ELISA) (Shanghai Jianglai Biotech, JL15386). The principle of ELISA kit was based on a double-antibody sandwich enzyme-linked immunosorbent assay to measure the concentration of ALDH18A1 by one step. Briefly, the plasma and the horseradish peroxidase (HRP)-labeled detection antibody were added to the microwells, which were pre-coated with ALDH18A1 antibody. Then, the mixture was incubate at 37°C for an hour. After washing with 1 × PBS, 3,3',5,5'-Tetramethylbenzidine was added as the chromogenic substrate. The wavelength at 450 nm excitation was used to assess the concentration of ALDH18A1 in the subject samples. The exact values were determined to the standard curve.

Data are presented as mean ± SEM. One-way ANOVA with Dunnett's multiple comparisons test was used and probability values P < 0.05 were considered statistically significant. The statistical analyses were performed with GraphPad Prism 5.0.

The patient II-2 was a 46-year-old woman in ongoing follow-up at the neurology outpatient clinic of the First Affiliated Hospital of Fujian Medical University due to a 20-year-long history of altered gait (Table 1). She was born in a non-consanguineous family and showed normal motor and psychological milestones. She experienced unstable walking and abnormal walking posture at the age of 26 and symptoms worsened after pregnancy. As gait abnormality progressed, she developed dysarthria, although no muscular atrophy or autonomic symptoms were present. No symptoms of cutis laxa were observed. Physical examination revealed decreased muscle strength and increased muscle tone of the lower limbs, hyperreflexia in all limbs, and bilaterally extensor plantar responses. Hoffmann sign was observed in her upper limbs. Spinal MRI presented atrophy of thoracic spinal cord (Supplementary Figure 1). Moreover, the plasma levels of amino acids, including proline, arginine, ornithine, and citrulline, were within the normal range (Supplementary Table 1). The proband had a younger brother, subject II-3, who carried the same mutations. He presented with difficulty in walking in his 30's. Walking abnormalities progressed gradually, but no other symptoms emerged. He can still walk independently at the age of 44. Their parents and elder sister carried the heterozygous mutation without detectable symptoms (Figure 1A).

We sequenced the exomes of the proband and her parents to 100X average depth in target regions, with 96% of the target base pairs on average having at least 20X coverage. WES analysis identified a heterozygous missense variant in ALDH18A1 (chr10: 95628421A>G hg38, NM_002860: exon8: c.880T>C: p.S294P), which was inherited from her mother (Figure 1B). This variant was completely absent from the ExAC, 1,000 Genomes Project, and Genome Aggregation Database (gnomAD) genome and exome databases. It also showed high evolutionary conservation (GERP score: 5.79), was predicted to be deleterious by at least three different algorithms (PolyPhen-2, MutationTaster, CADD), and was found to be located in the G5K domain of P5CS. Based on these results, we classified the variant as “Likely Pathogenic” according to ACMG/AMP guidelines (PM1, PM2, PM3, and PP1) (26).

As the clinical presentation of the proband was compatible with spastic paraplegia, and we did not find any other rare candidate variants in either the exon or conventional splicing region in the exome sequence, we then re-examined the exome data and found a potentially deleterious variant in the intron 1. This variant (chr10: 95653418T>C hg38, NM_002860: c.-28-13A>G) was intronic, completely absent from ExAC and 1,000 Genomes Project, occurred at extremely low frequencies among East Asian populations according to the gnomAD database (1/5,196), and has not been previously reported. This variant was inherited from her father and segregation in the family was compatible with pathogenicity of the variant. In silico analysis using SpliceAI, a novel tool for cryptic splice site prediction based on artificial neural network deep learning (15), predicted the disruption of an acceptor site located 13 bp downstream of the c.-28-13A>G variant. The loss of the original acceptor splice site likely led to the transcriptional skipping of exon 2.

We therefore sought to validate the SpliceAI prediction through analysis of cDNA reverse transcribed from RNA extracted from the peripheral blood samples of patients and family members. Two transcripts were observed in gel electrophoresis and Sanger sequencing confirmed the prediction that exon 2 was absent in the truncated transcript due to transcriptional skipping caused by the c.-28-13A>G variant (Figure 1C). The percentage of truncated transcript was 14%, which was quantified using HTS (Supplementary Table 2). This novel mRNA isoform lacked the original start codon necessary for initiation of full length ALDH18A1 transcription, instead initiating mRNA transcription at the second possible start codon, in a different downstream reading frame (Figure 1D). Each of these reading frames, i.e., for the correct and aberrant start codons, resulted in a different amino acid sequence, and hence one of which was non-functional P5CS protein.

The P5CS enzyme is known to catalyze proline and ornithine biosynthesis. The effects of disease-related mutations in ALDH18A1 are correlated with loss or decrease of P5CS function (2, 24, 27). To validate the functional involvement of ALDH18A1 in SPG9 in this family, we next performed ELISA analysis of the proband's plasma. The results showed that plasma P5CS concentration in patient II-2 was significantly decreased compared to that of the three, gender matched healthy controls (p < 0.05) (Figure 1E). Based on this functional validation, we proposed that compound heterozygosity for these two mutations most likely explains the phenotype of the patient.