Introduction

Most of the coding variants in the transthyretin gene (TTR; OMIM description: *176300; chromosome region 18q12.1) determine amyloidogenic processes, which cause TTR amyloidosis (OMIM phenotype description: #105210), an autosomal dominant hereditary lethal condition.1 The pathogenesis of this disorder is driven by the deposition of amyloid in extracellular spaces of several tissues and organs. This scattered distribution of TTR amyloid fibrils generally leads to a peripheral sensorimotor and autonomic neuropathy, exhibiting variability in cardiac, gastrointestinal, renal, and ocular involvement, as well as age of onset.2 The basis of this clinical variability may be partially explained by the high genetic heterogeneity, which has been linked to the existence of more than 100 amyloidogenic mutations in TTR gene.3 Some mutations are primarily associated with neuropathy, while others lead predominantly to cardiomyopathy, even though both clinical displays could be featured simultaneously to different extents.4 Moreover, although these symptoms may be present in patients carrying different TTR disease-causing mutations, clinical phenotypes are not always concordant, and the same mutation may be associated with different symptoms or signs, even within the same family.5 A clear example of the clinical heterogeneity of the carriers of the same TTR amyloidogenic mutation is the Val30Met mutation (rs28933979, c.148G>A, p.(Val50Met)), which is one of the leading causes of TTR amyloidosis.6 There are two European Val30Met endemic foci in Portugal and Sweden, and patients show two distinct clinical presentations. Val30Met Portuguese families experience early onset, strong severity, and high penetrance.7 Conversely, Val30Met Swedish patients have late onset, intermediate severity, and low penetrance.8 Moreover, high variability in penetrance according to the gender of the transmitting parent was observed in these geographic areas with a higher penetrance associated with maternal transmission of the Val30Met mutation.9 In literature there is also a case report regarding seven twin carriers of the Val30Met mutation with discordant clinical phenotypes, and the authors hypothesized genetic and epigenetic regulatory mechanisms as possible explanations.10 This strongly indicates that the basis of the complex genotype–phenotype correlation is not only due to the disease-causing mutation, but other genetic factors could have a role in the clinical presentation.11, 12, 13, 14 Our previous studies have suggested the putative impact of genetic variants located in the non-coding regions of the TTR gene in the modulation of the clinical phenotype of TTR amyloidosis. First, we observed the presence of genetic variants in binding sites for transcription factors involved in heart development in African populations, in accordance with the manifestation of cardiomyopathy in patients with African ancestry.15 Then, we analyzed the haplotypes associated with the Val122Ile (rs76992529, c.424G>A, p.(Val142Ile)) and Val30Met mutations, and identified independent haplotypes for the same mutation, suggesting multiple mutational events and different origins for these disease-causing mutations.16, 17 As non-coding variants located in regulatory elements can affect gene expression,18 we hypothesize that these variations influence the distribution and accumulation of TTR-derived fibrils across different tissues and organs, modulating consequently the disease phenotype.

Recently, the GTEx (Genotype-Tissue Expression) project (available at http://www.gtexportal.org/home/) has investigated the effect of coding and non-coding variants on tissue-specific gene expression.19 These effects of genetic variants can be used to estimate the genetically determined gene expression across human tissues, investigating its contribution with respect to disease pathogenesis.20 This approach permits to explore how genetic variation is associated with phenotypic presentation through its effect on gene expression in multiple tissues. In a recent study, we observed that genetically determined TTR expression varies across human populations and this variation appears to be in line with known epidemiological differences in the patients affected.21

On this basis, we evaluated the impact of the non-coding variants on the genetically determined TTR expression in Italian patients with TTR amyloidosis. Indeed, we hypothesize that non-coding variation regulating TTR gene expression in source and target tissues is one of the mechanisms involved in the phenotypic heterogeneity observed among the patients affected by TTR amyloidosis. However, we could not directly test the association of the non-coding variants with respect to every phenotypic presentation of the disease because we investigated a rare condition (Italian prevalence: 1 in 100 000 subjects) and the direct association analysis would have been underpowered. Conversely, we used an analytic approach that allowed us to combine the effect of non-coding variants on TTR expression across human tissues and to test the pattern of genetically predicted TTR expression with respect to disease symptoms.

Subjects and methods

Subjects

We recruited 129 symptomatic patients from six of the main Italian centers for the diagnosis and treatment of systemic amyloidosis: IRCCS-Istituto Auxologico Italiano – Piancavallo (VB), Careggi University Hospital – Florence, ‘San Giovanni Calibita’ Fatebenefratelli Hospital, Isola Tiberina – Rome, Fondazione Policlinico Universitario ‘A. Gemelli’ – Rome, ‘Vito Fazzi’ Hospital – Lecce, and Azienda Ospedaliera Universitaria Policlinico ‘G. Martino’ – Messina. The distribution of the samples across the recruiting centers is reported in Supplementary Table S1. This study was approved by the appropriate Bioethics Committees. TTR amyloidosis diagnosis was based on clinical signs/symptoms, amyloid deposits in biopsy samples, and the presence of an amyloidogenic TTR mutation. Clinical characteristics that describe the disease phenotype (ie, age of onset and organ involvements) were assessed as previously reported.22, 23, 24, 25, 26 Briefly, the following criteria were used to determine organ involvements: peripheral neurological involvement (peripheral neuropathy); cardiac involvement (electrocardiographic abnormalities); gastrointestinal involvement (gastric-emptying disorder and pseudo-obstruction); autonomic neurological involvement (orthostatic hypotension and voiding dysfunction); renal involvement (altered proteinuria and renal function); ocular involvement (keratoconjunctivitis sicca, secondary glaucoma, vitreous opacities, or pupillary abnormalities); and impotence (self-reported). Of the 129 patients, 55 agreed to provide the biological sample needed to conduct a re-sequencing analysis of their TTR gene, including coding and non-coding regions. As no criteria was applied to select the samples to be investigated in the sequencing analysis (except for the availability of the biological specimen needed for the DNA extraction), no selection bias is expected.

Genetic analysis

DNA, extracted through the standard phenol/chloroform technique, was checked by the quality control and its quantification. The quality control was carried out by spectrophotometry (Varian Cary 50 Bio UV/Visible Spectrophotometer, Varian Medical Systems, Palo Alto, CA, USA), while the quantification was carried out by Real Time-PCR (Quantifiler Duo DNA Quantification Kit, Applied Biosystems, Foster City, CA, USA). The TTR gene structure was analyzed through double-strand Sanger sequencing (3500 Genetic Analyzer). A 20 kb region (GRCh37 chromosome 18:29,164 000–29 185 000) was re-sequenced, including the TTR gene (7257 bp), ~6000 bp upstream (5′-UTR) and ~6000 bp downstream (3′-UTR) regions. Forty-five overlapping segments were amplified and sequenced using the same protocol reported in a previous study by Soares et al.27

Sequencing data analysis

We used Geneious v.9.1.5 software (http://www.geneious.com/ Biomatters Ltd., Auckland, New Zealand) to reconstruct the entire 20 kb region for each of the 55 patients. First, we checked the chromatograms by using the quality control and the trimming tools available through the software. Afterwards, we mapped the forward and reverse sequenced tracts with a reference sequence obtained from the 1000 Genomes database (GRCh37 chr18:29164000–29185000). To annotate the variants, we considered the following transcript and protein reference sequences: NM_000371.3 and NP_000362.1. The variant calling was performed with the Geneious v.9.1.5 software. The 53 identified variants (excluding the amyloidogenic mutations) were annotated by HaploReg v.4.1.28 The sequencing results were deposited in the Leiden Open Variation Database available at www.lovd.nl/TTR (Individual ID from 00103269 to 00103323).

Genetically determined TTR gene expression estimation

The GTEx Version 6 data (available at http://www.gtexportal.org/home/) were used as reference data sets for genetically determined gene expression. We extracted the comprehensive information (effect size and P-value) pertinent to 26 variants and 46 tissues (Supplementary Table S2). We then performed a linkage disequilibrium (LD) clumping analysis (P=1; r2=0.8; region size=10 kb) using Plink 1.0729 to obtain 46 tissue-specific data sets accounting for LD structure (Supplementary Table S3). All data sets consisted of genetic variants in non-coding regions of the TTR gene, with the exception of the missense not-pathological rs1800458. We used these tissue-specific clumped data sets to calculate the individual polygenic scores for the genetically determined TTR expression, on the basis of effect allele count and allele effect size, using Plink 1.07. Specifically, for each tissue the scores were calculated weighting the allele dosage by its effect on TTR expression in the corresponding tissue.

Statistical analysis

The initial analysis was conducted on the phenotypic information available for the 129 symptomatic patients. A hierarchical clustering analysis, using the Jaccard method for binary variables, was conducted to correlate the amyloidogenic mutation and the clinical information. This analysis was performed using ‘arules’ R package30 and the cluster significance using ‘pvclust’ R package (available at https://cran.r-project.org/web/packages/pvclust/index.html) by applying 100 000 bootstrap replications. We considered clusters with an approximately unbiased P-value >95% as significant.

The second phase of the study was to investigate the relationship between coding and non-coding TTR variation and clinical characteristics of the 55 sequenced patients. The haplotype reconstruction (for all the identified variants and amyloidogenic mutations) was conducted using PHASE 2.131 and the haplotype association with clinical phenotypes was performed with Haploview 4.2.32 We performed a hierarchical clustering analysis on the polygenic scores using MeV 4.8.1 software (available at http://www.tm4.org/mev.html) and calculated the cluster significance using ‘pvclust’ R package as described above. We generated a distance metric using the Spearman rank correlation method, and the clustering was conducted using the average linkage method. For the significant clusters, we checked the statistical significance of the clinical phenotypes by comparing their observed frequencies with the null distribution generated through 10 000 random permutations of the phenotypic data with respect to polygenic scores.

Results

The characteristics related to the amyloidogenic mutations and clinical phenotypes of the 129 Italian patients are shown in Table 1. The most frequent mutation is Val30Met (N=55), followed by Phe64Leu (rs121918091, c.250T>C, p.(Phe84Leu); N=28), Glu89Gln (rs121918082, c.325G>C, p.(Glu109Gln); N=28), Thr49Ala (rs121918081, c.205A>G, p.(Thr69Ala); N=7), Ala120Ser (rs876658108, c.418G>T, p.(Ala140Ser); N=3), Ile68Leu (rs121918085, c.262A>T, p.(Ile88Leu); N=3), Val122Ile (N=3), Ala109Ser (rs267607159, c.385G>A, p.(Ala129Thr); N=1), and Thr59Lys (rs730881163, c.236C>A, p.(Thr79Lys); N=1). We observed similar symptom patterns across patients with different amyloidogenic mutations with few exceptions (Figure 1). Cardiac, autonomic neurological, peripheral neurological, and gastrointestinal symptoms were observed in most of the amyloidogenic mutations. However, we observed few strong correlations between amyloidogenic mutations and some phenotypic traits. Ocular symptoms were only observed in Val30Met patients. Ile68Leu and Val122Ile patients almost exclusively showed the cardiac phenotype. Applying hierarchical clustering (Figure 2), we observed significant correlation between the Val30Met mutation and ocular symptoms (approximately unbiased P>95%) and between the Phe64Leu mutation and renal involvement (approximately unbiased P>95%). Both the Ile68Leu and Val122Ile mutations were significantly associated with the cardiac symptoms (approximately unbiased P>95%).

Table 1 Description of the clinical features and the amyloidogenic mutations of the 129 Italian patients
Full size table
Figure 1
figure 1

Radar diagram of symptoms across amyloidogenic mutations. The percentages are calculated with respect to the total number of organ involvements reported for each TTR mutation. The full colour version of this figure is available at European Journal of Human Genetics online.

Full size image
Figure 2
figure 2

Hierarchical clustering based on clinical symptoms and amyloidogenic mutation reported in the 129 patients with TTR amyloidosis. Red boxes highlight the significant clusters (approximately unbiased (AU) P>95%; bootstrap probability (BP)). Gastr.Int, gastrointestinal involvement; SNA, autonomic neurological involvement; SNP, peripheral neurological involvement. The full colour version of this figure is available at European Journal of Human Genetics online.

Full size image

Re-sequencing coding and non-coding regions of the TTR gene in 55 patients, we identified 53 variants (Supplementary Table S4), including 43 SNPs and 10 INDELs (21 in the upstream region, 8 in introns, 1 in the coding regions, and 24 in the downstream region). Four of the 53 variants (NC_000018.9:g.29164968C>T, g.29179948G>A, g.29182334A>G, and g.29182576T>C) are novel and not annotated in the international reference databases and genome browsers (dbSNP, 1000 Genomes, UCSC Genome Browser, Ensembl). Functional annotation analysis (Figure 3) indicated that 9 variants are associated with regulatory elements in the promoter region, 31 are associated with regulatory elements in the enhancer regions, 20 are found in DNAse-sensitive loci, 5 are located in protein-binding sites, 39 alter transcription factor-binding sites, and 5 are known expression quantitative trait loci.

Figure 3
figure 3

Functional annotation of the variants identified by TTR gene re-sequencing.

Full size image

We performed the haplotype inference considering the 53 identified variants and the 7 amyloidogenic mutations. The haplotypes associated with the disease-causing mutations of the 55 patients are reported in Supplementary Table S5. Through the haplotype association analysis, we found a specific Val30Met haplotype significantly associated with ocular symptoms (P=8 × 10−4). This Val30Met haplotype contains two regulatory variants: rs5823817 and rs58272364. These variants are intronic and are associated with regulatory elements in enhancer regions. Furthermore, rs5823817 is also located in in a DNAse-sensitive locus.

As explained in Subjects and Methods, we used our sequencing data to estimate genetically determined TTR expression across human tissues. We applied a hierarchical clustering analysis to these data to identify subjects with similar TTR expression (Figure 4). We identified five significant patient clusters (approximately unbiased P>95%), which include patients with different amyloidogenic mutations: cluster 1 (Val30Met, Phe64Leu, and Ile68Leu); cluster 2 (Val30Met and Glu89Gln); cluster 3 (Phe64Leu and Val122Ile); cluster 4 (Ala120Ser and Phe64Leu); and cluster 5 (Val30Met and Phe64Leu). We performed a permutation analysis to verify whether the patients included in the clusters identified on the basis of the genetically predicted TTR expression are associated with specific symptoms. Clusters 3 and 4 were excluded from further analysis because of their small sample size (N=3). Subjects included in cluster 1 showed a higher age of onset and a lower prevalence of autonomic neurological symptoms than those expected by chance (ppermutation=0.013 and 0.035, respectively; Supplementary Figure S1). An increased prevalence of autonomic neurological symptoms (trend significance; ppermutation=0.087) was observed in the subjects included in cluster 2 (Supplementary Figure S2). In cluster 5, we observed subjects with a higher prevalence of gastrointestinal phenotype than that expected by chance (ppermutation=0.047; Supplementary Figure S3).

Figure 4
figure 4

Hierarchical clustering using tissue-specific polygenic score of the 55 patients re-sequenced. On the top the 55 patients and on the left the 46 considered tissues. Red boxes highlight the significant clusters (approximately unbiased P>95%). Red, green, and black colors correspond to genetically determined overexpression, under-expression, and little differential expression, respectively. The color scheme is also reported at the top of the figure. The full colour version of this figure is available at European Journal of Human Genetics online.

Full size image

Discussion

Our findings provided novel insights into the molecular mechanisms at the basis of the genotype–phenotype correlation of TTR amyloidosis. Considering the clinical information available for the 129 Italian patients investigated in the present study, we confirmed a strong phenotypic heterogeneity: the main symptoms (neurological, cardiac, and gastrointestinal) were observed across different amyloidogenic mutations. This indicates that the symptom pattern in TTR amyloidosis is not exclusively determined by the disease-causing mutation but other mechanisms should be involved. However, we also observed few strong genotype–phenotype correlations: Val30Met mutation with the ocular involvement; and Val122Ile and Ile68Leu mutations with the cardiac phenotype. These phenotype–genotype correlations are supported by previous studies. The ocular involvement in patients with TTR amyloidosis is generally a rare manifestation33 because it depends on numerous factors (eg, liver transplantation, disease age of onset, and disease duration), but it is well described in several studies of Val30Met TTR amyloidosis in Portuguese families.34, 35 The Val122Ile mutation is widely described in the literature as the leading cause of cardiac amyloidosis in patients with African ancestry.36 This correlation was also confirmed in patients of European descent.25, 37 The Ile68Leu mutation was previously described in Italian patients with cardiac or mixed neurological/cardiac phenotypes.4, 25

The re-sequencing analysis of the TTR gene (coding and non-coding regions) in Italian patients allowed us to deepen our previous finding on the role of non-coding TTR variations in the disease genotype–phenotype correlation.15, 16, 17 Investigating the TTR haplotypic structure, we observed that the correlation between the Val30Met mutation and the ocular manifestation can be attributed to a specific Val30Met haplotype, which includes two regulatory non-coding variants. As also reported in genetic studies of different phenotypic traits,38, 39 the ocular manifestation in TTR amylodosis in our sample is due to a combination of the effects of a coding amylodogenic mutation and SNPs located in non-coding elements.

Using of the data from the GTEx consortium,19 we leveraged our sequencing data to calculate the genetically predicted expression profiles of the TTR gene across 46 human tissues. These expression profiles reflect the patterns of TTR expression influenced by the genetic control. Consequently, they should be related to potential pathogenic mechanisms linked to TTR expression across human tissues and not to the consequences of the disease progression. Considering these expression profiles, we identified three patient clusters with specific expression patterns that were associated with particular phenotypic presentations: cluster 1, late onset and low prevalence of autonomic neurological involvement; cluster 2, high prevalence of autonomic neurological involvement; and cluster 5, high prevalence of gastrointestinal symptoms. This indicates that the genetic variation involved in the regulation of TTR expression across human tissues is one of the mechanisms responsible for the clinical presentations observed in the carriers of TTR amyloidogenic mutations. Furthermore, the effects of non-coding regulatory elements appear independent of the amyloidogenic mutations: variants involved in TTR gene expression regulation are involved in the disease manifestation with similar effects across different disease-causing mutations.

Owing to the limited availability of tissue samples taken from affected patients, mechanisms related to gene expression in TTR amyloidosis have been poorly investigated. Recently, an in vitro study demonstrated that Schwann cells can contribute to neurodegeneration through the local expression of mutated TTR.40 TTR gene expression patterns across different tissues, including source and target organs, could contribute to the symptoms observed in patients, which is in agreement with our current findings.

In conclusion, the present study provided novel data regarding the role of non-coding variation and the expression profiles of TTR gene in affected patients, also introducing an innovative approach to investigate the mechanisms at the basis of disease symptoms. Further experimental studies are needed to confirm our results. Furthermore, other mechanisms (eg, epigenetics and modifier genes) can also contribute to the disease genotype–phenotype correlation, and local gene expression profile may reflect one of the molecular processes that contribute to the genotype–phenotype correlation of TTR amyloidosis.