This article is mentioned in:

Introduction

Next-generation sequencing (NGS) technologies with massively parallel sequencing are widely used for medical genomics studies (1). NGS techniques, such as whole-exome sequencing (WES) and whole-genome sequencing (WGS), are more desirable than individual gene sequencing due to high coverage sequencing at a lower cost (2). However, data processing and analysis remain a limitation in WES and WGS (1). Specifically, identifying all single nucleotide variants (SNVs) and short insertion/deletions (indels) from the protein-coding region is a challenge with WES data analysis (3). Confounding factors such as DNA quality and numerous potential errors during the library preparation, DNA sequencing, alignment and mapping steps affect the accuracy of the variants called (3). Multiple quality control steps are employed in a standard bioinformatics pipeline; however, false positive and negative variants still occur (4). Therefore, visual inspection of read alignments is key to accurately identify the variants from NGS data (4,5).

Autosomal dominant polycystic kidney disease (ADPKD) is the most common type of inherited cystic kidney disease, with an estimated prevalence of 9.3 per 10,000 people worldwide (6) and is characterized by development of multiple cysts in both kidneys. Due to the enlargement of the kidneys and progressive loss of renal function, ~50% of patients with ADPKD suffer from end-stage renal disease (ESRD) by age 60 (7). ADPKD is primarily caused by mutations in the polycystin 1, transient receptor potential channel interacting (PKD)1 and PKD2 genes (7). PKD1 is composed of 46 exons with a coding length of 12,912 bp (NM_001009944.3) (8), whereas PKD2 is composed of 15 exons with a coding length of 2,907 bp (NM_000297.4) (9). Besides being a large gene, sequencing PKD1 is complicated by the presence of six pseudogenes (PKD1P1-PKD1P6) that share >97% sequence similarity with exons 1–33 of PKD1 (10).

Molecular analysis of ADPKD is performed using several techniques such as long-range PCR followed by direct Sanger sequencing (11), multiple ligation probe assay (12), NGS techniques (13,14) or a mixture of the aforementioned approaches (15). Conventionally, long-range PCR is used to exclude pseudogenes (11). However, the subsequent Sanger sequencing for PKD1 is laborious, expensive and time-consuming (16). Due to these factors, NGS technologies such as WGS (13) and WES (14), have increasingly been utilized to genotype PKD1 and PKD2 in patients with ADPKD. Compared with Sanger sequencing, a recent study reported that WGS has 100% sensitivity and specificity in detecting variants associated with ADPKD (13). However, the utility of WES for ADPKD remains in unknown, as recent study reported the sensitivity to be limited at 50% (14).

Our previous studies reported successful use of WES to identify genetic mutations for several types of monogenic disease that belong to the group of inborn errors of immunity (17,18). In the present study, the application of WES to ADPKD was evaluated. Key visual inspection steps in identifying a novel insertion mutation in PKD1, which was not identified by the variant caller, are described.

Subjects and methods

Study subject

A 50-year-old woman was diagnosed with advanced chronic kidney disease in September 2013 during health screening in Kuala Lumpur Hospital, Malaysia. Ultrasound was performed because the patient had abnormal kidney function. Ultrasonography of her kidneys demonstrated bilateral polycystic kidneys with features which were highly suggestive of ADPKD. Diagnosis of ADPKD was made clinically based on the ultrasound result in September 2013. The patient was then recruited into the research cohort before her kidney transplant in August 2015.

WES

Genomic DNA was extracted from whole blood in EDTA tubes using QIAsymphony DSP DNA Midi kit (cat. no. 937255; Qiagen GmbH) on a QIAsymphony SP instrument (Qiagen GmbH). The DNA concentration and purity were evaluated through optical density measurement at 260 nm and 260/280 ratio respectively, using the QIAxpert Slide-40 (cat. no. 990700; Qiagen GmbH) on a QIAxpert System (Qiagen, GmbH). Two µg of genomic DNA was fragmented into 150–200 bp using Covaris LE220-plus focused-ultrasonicator (Covaris, Inc.). A genomic library with fragment size of ~330 bp was constructed using SureSelectXT Reagent kit (cat. no. G9641C; Agilent Technologies, Inc.). Exome enrichment was performed using the SureSelect Human All Exon V6 kit (cat. no. 5190-8864; Agilent Technologies, Inc.) with a target size of 60 Mb. The size of PCR enriched fragments was verified using Agilent DNA 1000 kit (cat. no. 5067-1504; Agilent Technologies, Inc.) on a 2100 Bioanalyzer instrument (Agilent Technologies, Inc.). The final library was quantified using qPCR according to the Illumina qPCR Quantification Protocol Guide (KAPA SYBR FAST qPCR Master Mix (2X) Universal; cat. no. KK4602; Kapa Biosystems, Inc.). The loading concentration of the final library was 300 pM. Paired-end reads of 2×101 bp were sequenced using the HiSeq 3000/4000 SBS kit (300 cycles) (cat. no. FC-410-1003; Illumina, Inc.) on a HiSeq 4000 System (Illumina, Inc.) with a minimum coverage of 100×. The raw data were converted into FASTQ format.

Bioinformatics analysis

The bioinformatics processing pipeline for germline short variant discovery was modified from the Genome Analysis Toolkit (GATK) Best Practices Workflows (version 4.1.2.0) (19). Briefly, pre-processing of the FASTQ file began with addition of specific read group information and tagging of Illumina adapters using Picard (version 2.20.1) (20). Next, the reads were aligned and mapped to the human reference genome GRCh38 (GCA_000001405.15) using the Burrows-Wheeler Aligner-maximal exact matches (BWA-MEM) (version 0.7.17-r1188) (21). Additional quality control step in BAM file processing was included to unmap contaminant reads using Picard (version 2.20.1) (20). Following alignment, duplicate reads were marked using Picard (version 2.20.1) (20) and base quality score recalibration (BQSR) was performed using GATK (version 4.1.2.0) (19). Variants including SNVs and indels were called using HaplotypeCaller (version 4.1.2.0) (22). Finally, the resulting variant call format file was annotated using the web-based wANNOVAR tool (accessed in March 2020) (23).

The alignment of reads in the binary alignment and map (BAM) file was visualized using the Integrative Genomics Viewer (IGV) (version 2.8.10) from The Eli and Edythe L. Broad Institute of MIT and Harvard (24). Sequence similarity analysis was performed using the Basic Local Alignment Search Tool (BLAST) (25). The pathogenicity of variants was computationally evaluated using in silico prediction tools, namely Sorting Intolerant from Tolerant (version 2.3) (26), MutationTaster (version 2) (27) and PolyPhen-2 (version 2.2.2) (28). Nucleotide and protein changes were identified by comparison with National Center for Biotechnology Information (NCBI) reference sequences of PKD1 (NM_001009944.3) (8) and PKD2 (NM_000297.4) (9). The detected mutation sites were compared with the ADPKD Variant Database (PKDB) (29) and NCBI Single Nucleotide Polymorphism Database (dbSNP) human build 155 (30). Allelic frequency of the variants was checked using the genome Aggregation Database (gnomAD) (31). The final decision on pathogenicity of the detected mutations was based on the American College of Medical Genetics and Genomics (ACMG) classification (32). Changes to protein sequences following an indel event were determined using the ExPASy translation tool (33).

Mutation validation

Primers were designed to flank the targeted region using Primer3 (version 0.4.0) (34) and validated using Primer-BLAST (35). Using genomic DNA extracted from whole blood of patient, primer sequences (forward, 5′-CTGCTCTTCCTGCTTTTGGT-3′ and reverse, 5′-CCGTACCCACCTCCTTGAC-3′) were used to amplify a product of 633 bp from the PKD1 gene using MyFi™ Mix kit (Meridian Bioscience, Inc.). Genomic DNA from a healthy unrelated individual was used as the control. The PCR cycling conditions were initial denaturation at 95°C for 3 min, followed by 33 cycles of denaturation at 95°C for 30 sec, annealing at 60°C for 30 sec and extension at 72°C for 40 sec and final extension at 72°C for 8 min. PCR products were purified and subjected to bidirectional Sanger sequencing using BigDye™ Terminator v3.1 Cycle Sequencing kit (Applied Biosystems; Thermo Fisher Scientific, Inc.) on a 3500×l Genetic Analyzer (Applied Biosystems; Thermo Fisher Scientific, Inc.). The chromatograms were viewed and analyzed using FinchTV (version 1.4.0; Geospiza, Inc.).

Results

Case presentation

The patient was the sixth of eight children. Her father had ischemic heart disease and her mother had diabetes mellitus. No family history of renal disease was noted; however, only one of her brothers was screened for the disease. The patient developed ESRD in September 2013 and hemodialysis was initiated. Her kidney ultrasound demonstrated >15 cysts in each kidney with bipolar lengths of 15 and 17 cm. Diagnosis of ADPKD was made despite negative family history and genetic testing. Bilateral nephrectomy was performed in November 2014 and she subsequently received a living kidney transplant from her younger brother in August 2015. The kidney transplant was successful with good kidney function.

Bioinformatics analysis of WES

WES generated 70,808,986 paired-end reads, of which 69,941,464 reads (98.77%) were properly paired and mapped to the human reference genome GRCh38. A total of 110,294 reads suspected of cross-species contamination due to extremely short alignments with clipping on both sides were unmapped. Variant calling annotated 24,901 SNVs and 703 indels in the exonic and splice site regions. A total of 13 SNVs was identified in the PKD1 and PKD2 genes. Of these, eight were synonymous and five were non-synonymous mutations. The pathogenicity and allelic frequency for all variants were carefully examined. All five non-synonymous mutations were assessed to be benign using in silico prediction tools and were present in >5% of the population. The workflow for bioinformatics processing and variant prioritization is presented in Fig. 1.

Figure 1. Workflow of WES for the patient. WES is divided into bioinformatics processing and variant prioritization. In bioinformatics processing, a total of 70,808,986 paired-end reads were generated, of which 69,941,464 reads were properly paired and mapped to the human reference genome. Subsequently, 110,294 reads suspected as cross-species contamination were unmapped as part of quality control. The final BAM file was subjected to variant calling and annotated 24,901 SNVs and 703 indels. In the variant prioritization step, variants were filtered against PKD1 and PKD2 genes. A total of 13 SNVs was identified, of which eight were synonymous; five were non-synonymous mutations. The pathogenicity of these five SNVs was predicted using in silico prediction tools SIFT, MutationTaster and PolyPhen-2. All five SNVs were assessed to be benign. BAM, binary alignment and map; SNV, single nucleotide variant; indel, insertions/deletion; SIFT, Sorting Intolerant from Tolerant.

Identification of atypically aligned soft-clipped reads

As no pathogenic genetic variants were demonstrated, the read alignments of PKD1 and PKD2 genes were visually inspected using IGV. One position at chr16:2,090,586 of exon 45 of PKD1 was assessed as having numerous reads containing long sequences of soft-clipped bases (Fig. 2). This region was unique because all soft-clipped bases were in alignment with each other, which prompted further evaluation. This atypical alignment may have resulted from strand bias during capture, poor mapping, sequencing errors or DNA template contamination.

Figure 2. IGV of atypically aligned soft-clipped reads. (A) IGV visualization of short reads alignment using the binary alignment and map file. Pink, forward read; blue, reverse read. Mismatches and soft-clipped bases are bases that do not match the human reference genome and were presented as colored letters on the read. No letters were presented if the reads fully matched the human reference genome. Upon visual inspection, one region at the beginning of exon 45 in PKD1 contained numerous reads with long soft-clipped bases. Red arrow indicates location of a potential insertion mutation. (B) Single base mismatch. In this example, certain reads contained base G, which was matched to the human reference genome. However, certain reads contained base A, which was a mismatch. This was an example of a heterozygous G/A single nucleotide variant detection. (C) Section of ten reads with zero mismatches. No letters were presented on reads as they fully matched the human reference genome. White space indicates interval between aligned reads. (D) Section of five reads with soft-clipped bases. The colored letters indicate this section of the reads was mismatched with the human reference genome. Long soft-clipped bases were present in both forward and reverse reads and mismatches aligned with one another. IGV, Integrative Genomics Viewer.

Soft-clipped reads are of high-quality mapping and good base score

The number of reads that spanned the affected site was counted (Fig. 3A); there was a total of 91 reads, with 52 reads without any soft-clipped bases and 39 reads with soft-clipped bases. The forward reads were colored pink and the reverse reads blue. For 52 reads without soft-clipped bases, 25 were forward and 27 were reverse reads. For the 39 reads with soft-clipped bases, 21 were forward and 18 were reverse reads. The coverage of reads at the site of interest was sufficiently high to suggest heterozygosity. There was no evidence of strand bias as the distribution of forward and reverse reads was balanced.

Figure 3. Flow chart of assessment of atypical alignment of soft-clipped reads using IGV and BLAST. (A) Determination of quantity and coverage of reads with and without soft-clipped bases. The total coverage of 91 was high and equally distributed between the with and without soft-clipped groups, which suggested true heterozygosity. (B) Determination of MAPQ of reads and QV of bases. MAPQ=60 indicated high confidence that the read was correctly aligned; MAPQ=0 indicated the read could be mapped to more than one location in the human reference genome. The potential QV scores for bases ranged from 6 to 43. Higher score indicated higher probability of correctly calling the base. In the example presented, the read with soft-clipped bases had MAPQ=60 (red box) and the soft-clipped base G had QV 40 (purple box). (C) Checking mapping quality of mates of soft-clipped reads. If mates had low-quality MAPQ=0, this indicated potential sequencing artifacts or errors during library preparation. In the example presented, the mate (left) for the read presented (right) had MAPQ=60 (green box). (D) Determination of origin of the soft-clipped bases using IGV built-in BLAT function and BLAST. Atypical alignment with soft-clipped bases results from non-human DNA contamination, which can be determined from BLAST results. In the example presented, the results obtained from BLAST demonstrated that soft-clipped read had >97% identity (orange box) with Homo sapiens polycystin 1, which indicated no evidence of cross-species contamination. When all conditions were fulfilled, the atypical alignment likely represented the actual genomic sequence. IGV, Integrative Genomics Viewer; BLAST, Basic Local Alignment Search Tool; BLAT, BLAST-like Alignment Tool; MAPQ, mapping quality; QV, quality value; NCBI, National Center for Biotechnology Information.

The mapping quality (MAPQ) of the soft-clipped reads and the quality value (QV) of the soft-clipped bases were checked. In the example presented in Fig. 3B, the read was 101 bp in length with 63 bases matched and 38 bases soft-clipped on the right. After BWA-MEM, the MAPQ of reads ranged from 0 to 60. This read had the highest possible value of MAPQ=60, which indicated high confidence that the read was correctly aligned. Read with low-quality mapping (MAPQ=0) appeared translucent on IGV and indicated that the read could be mapped to more than one location in the human reference genome. Assessment of all 91 reads that spanned the affected region demonstrated that all reads had MAPQ=60. Following BQSR, the QV ranged from a minimum of 6 to a maximum of 43. The base G in the example presented had a QV score of 40. Generally, soft-clipped bases had high QV scores of 30–42. These findings ruled out poor mapping and sequencing errors as potential causes for atypical alignment.

Subsequently, the mapping parameters of mates of the soft-clipped reads were checked using the IGV option to view the reads as pairs. The mate of the read presented in Fig. 3B was used as the example presented in Fig. 3C. The right read had 38 bases soft-clipped and the left read was 101 bases and fully matched with MAPQ=60. Mates with MAPQ=0 may have been due to sequencing artifacts or errors during the library preparation. While in ‘view as pairs’ mode, the start and end of insert size were compared. Read pairs with the same start and end positions may indicate duplication events that were not removed during quality control and may have given a false sense of coverage. No evidence of duplication events was discovered and all mates of soft-clipped reads had MAPQ=60.

Reads with soft-clipped bases are of PKD1 origin

The final step was to identify the origin of soft-clipped bases. Atypical alignment with soft-clipped bases may result from non-human DNA contamination (36). To evaluate this, BLAST was performed on five reads with long soft-clipped bases (Fig. 3D); reads had >97% similarity to PKD1. This was also performed using IGV built-in BLAST-like Alignment Tool (BLAT) function on the soft-clipped bases to assess whether the soft-clipped sequence could be mapped to other parts of the human reference genome. Five results from BLAT analysis demonstrated that the soft-clipped bases were partially matched to PKD1. The results demonstrated the atypically aligned soft-clipped reads were of PKD1 origin.

Sanger sequencing confirms presence of a 20-bp novel insertion within the soft-clipped region

The atypical alignment at PKD1 was hypothesized to be caused by an indel mutation. Sanger sequencing confirmed the presence of two SNVs (rs1384786564 and rs1337808849) and a 20-bp insertion (c.12143_12144insTCCCCGCAGTCTTCCCCGCA; Fig. 4A). This 20-bp insertion was considered novel because it has not been previously reported in commonly used databases including PKDB, dbSNP (build 155) and gnomAD. This mutation was predicted to change the reading frame from codon 4,048 for 156 amino acids, followed by introduction of a premature stop codon at codon 4,204. The prematurely truncated protein was predicted to be 100 amino acids shorter than the wild-type protein of 4,303 amino acids (Fig. 4B). Assessment using ACMG classified this novel 20-bp insertion as pathogenic.

Figure 4. Identification of a novel PKD1 insertion in the patient. (A) A total of two SNVs and a novel 20-bp insertion, leading to a frameshift and premature stop codon 157 amino acids downstream, were confirmed in the patient using Sanger sequencing. These mutations were absent in the control. Three red arrows indicated the mutation sites. Red arrow 1 indicates a synonymous mutation c.12141C>T; p.Leu4047= (rs1384786564). Red arrow 2 indicates a non-synonymous mutation c.12142G>C; p.Val4048Leu (rs1337808849). Red arrow 3 indicates a novel 20-bp insertion c.12143_12144insTCCCCGCAGTCTTCCCCGCA; p.Val4048LeufsTer157. (B) c.12143_12144insTCCCCGCAGTCTTCCCCGCA mutation and rs1337808849 in exon 45 of PKD1 were predicted to change the reading frame and introduced an early stop codon (TGA) at position 4204 (p.Val4048LeufsTer157). This would lead to a prematurely truncated protein 100 amino acids shorter than the wild-type protein of 4,303 amino acids. Red arrow 4 indicates location of the 20-bp insertion into the wild-type sequence. SNV, single nucleotide variant; CDS, coding sequence.

Discussion

A novel 20-bp insertion mutation in PKD1 from a patient with ADPKD was identified using WES. However, this mutation was not identified by the variant caller. Instead, a region with atypical alignment of numerous soft-clipped reads was identified by visual inspection of read alignments. A total of four visual inspection steps was outlined using IGV and BLAST to rule out potential errors that may yield soft-clipped reads. Finally, the hypothesis that an indel event caused the atypically aligned soft-clipped reads was confirmed by Sanger sequencing. To the best of our knowledge, this 20-bp insertion in PKD1 is novel and predicted to introduce a premature stop codon, acting as a loss-of-function mutation. Despite the lack of ADPKD in the family history, assessment using ACMG classified this mutation as pathogenic, which supported the clinical diagnosis of ADPKD for the patient.

SNV and indel detection using WES is routine; however, researchers rely on the output generated by automated pipelines and inaccurate or false variant calls occur. Visual inspection of aligned reads is a key step for variant discovery and validation (4,5). The proposed visual inspection steps can be incorporated into WES data analysis workflow following variant prioritization (Fig. 5). Visual inspection at the genes of interest may reveal atypical alignment. It has been reported that atypical alignment with soft-clipped bases results from contamination with non-human DNA (36). Samson et al (36) reported that non-human DNA originating from the oral microbiome aligns with the human reference genome, producing atypical alignment and false positive variants. However, the DNA source of the aforementioned study was saliva, whereas the DNA source in the present study was peripheral blood. Visual inspection using BLAST and BLAT confirmed that the atypical alignment in the present patient was not due to cross-species contamination.

Figure 5. Proposed incorporation of visual inspection step in WES data analysis workflow. Typical WES data analysis workflow begins with alignment and mapping of short reads to the human reference genome. This is followed by variant calling and annotation. The annotated variants are filtered and prioritized according to the disease of interest and associated genes. The final step is variant validation. Visual inspection may be incorporated following variant prioritization. *, step in which the 20-bp insertion was missed by the variant caller. WES, whole-exome sequencing.

In the present study, an indel mutation leading to truncation in PKD1 was found in a patient with ADPKD. Generally, mutations are more commonly identified in PKD1 (~85%) than PKD2 (~15%) (37–39). At the time of analysis, a total of 1,225 PKD1 and 196 PKD2 pathogenic variants were cataloged in the PKDB (29). For PKD1, indels that lead to a truncated protein are the most common mutation (36.7%; 449/1,225), followed by missense (22.6%; 277/1,225), nonsense (21.5%; 263/1,225) and other forms of mutation (19.3%; 236/1,225).

Despite numerous variants reported in PKD1 and PKD2, no pathogenic variant has been reported in these two genes for 11–39% of patients with ADPKD in large cohort studies (37–39). For these patients, several other genes have been reported to cause ADPKD, such as glucosidase II alpha subunit (40), DnaJ heat shock protein family 40 member B11 (DNAJB11) (41), ALG8 alpha-1,3-glucosyltransferase and protein kinase C substrate 80K-H (39). Furthermore, numerous other candidate genes have been reported to affect progression and severity of ADPKD. Using WES, Hu et al (42) assessed 313 genes associated with polycystic kidney disease; demonstrated that molecular analysis of patients with ADPKD using global approaches such as WES, was more useful than individual gene sequencing using Sanger sequencing. For patients without mutations in PKD1 and PKD2, WES data can be used to screen for other cystic genes that may point to atypical ADPKD (16).

Sequencing PKD1 is challenging due to the presence of six pseudogenes (10). However, WGS is reported to overcome this (43) and has 100% sensitivity in detecting pathogenic variants in PKD1 (13). To the best of our knowledge, however, the diagnostic ability of WES for PKD1 is still uncertain. Al-Muhanna et al (44) reported that 100% coverage of PKD1 is possible with WES, whereas Ali et al (14) reported poor coverage of PKD1 at duplicated regions. The difference in reported sensitivity of WES may be due to different capture kits and data analysis methods. NGS is becoming more affordable; however, WGS is considerably more expensive compared with WES (2). WES that specifically targets protein-coding exomes at 100× coverage generates ~6 GB of data, whereas WGS that targets the entire genome at 30× coverage generates ~90 GB of data (1). The markedly larger data output from WGS demands higher computing power for data analysis, which is not readily available in many laboratories. Hence, WES is a cheaper and faster solution than WGS.

In clinical practice, diagnosis of ADPKD is usually made using ultrasonography, together with clinical presentations and family history. However, there is lower sensitivity of diagnosis using ultrasonography for younger individuals (45). Therefore, genetic testing may be beneficial to provide a definite ADPKD diagnosis. Genetic testing may also be used for testing patients with negative family history and atypical clinical presentation and for the selection of family members for living kidney transplantation (46).

In the present study, a novel 20-bp insertion in exon 45 of PKD1 in a patient with ADPKD was demonstrated using WES. This frameshift mutation was predicted to be pathogenic due to the introduction of a premature stop codon. Therefore, WES was demonstrated to be a suitable option for genetic testing for ADPKD, which is cheaper and faster than WGS or Sanger sequencing. However, data analysis of WES relying on a standard bioinformatics pipeline may miss out disease-causing variants. As demonstrated in the present study, visual inspection of read alignments was key to identifying the pathogenic variant. The proposed visual inspection steps can be incorporated into a typical WES data analysis workflow and may improve diagnostic yield. The limitations of this study include the requirement of experienced personnel to perform visual inspection in WES data analysis, and the results were based on a single case. Therefore, future research is warranted to validate the use of visual inspection on a larger patient cohort.

Acknowledgements

The authors would like to thank the Director General of Health, Malaysia for permission to publish this article. The authors would also thank Mr Muhammad Johari (Allergy and Immunology Research Centre, Institute for Medical Research) for assistance in DNA extraction.

Funding

The present study was funded by a grant from the Ministry of Health, Malaysia (grant no. NMRR-17-892-35929).

Availability of data and materials

The datasets generated and/or analyzed during the current study are available in the NCBI Sequence Read Archive under accession number PRJNA798582 (ncbi.nlm.nih.gov/sra/?term=PRJNA798582). The human reference genome GRCh38 (filename: GCA_000001405.15_GRCh38_no_alt_plus_hs38d1_analysis_set) can be downloaded from https://ftp.ncbi.nlm.nih.gov/genomes/all/GCF/000/001/405/GCF_000001405.26_GRCh38/GRCh38_major_release_seqs_for_alignment_pipelines/. dbSNP (build 155) (filename: GCF_000001405.39.gz) can be downloaded from https://ftp.ncbi.nlm.nih.gov/snp/archive/b155/VCF/.

Author's contributions

BTK conceived the study, acquired the funding, performed laboratory experiments and data analysis, interpreted the data and drafted the manuscript. MYC performed laboratory experiments and drafted the manuscript. JI and NKF recruited the patient, interpreted the data and constructed figures. SYY provided the clinical history of the patient and interpreted the data. NM and MA conceived the study and acquired the funding. AMR and SBM designed and supervised the study. BTK, AMR and SBM confirm the authenticity of all raw data. All authors have read and approved the final manuscript.

Ethics approval and consent to participate

Written informed consent for participation was collected from the patient prior to blood sample collection. The present study was approved by the Medical Research and Ethics Committee, Ministry of Health, Malaysia (approval no. KKM/NIHSEC/P17-1084) and was performed according to the Declaration of Helsinki.

Patient consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Glossary

Abbreviations

Abbreviations:

References