Heterozygous Tropomodulin 3 mice have improved lung vascularization after chronic hypoxia

Abstract

The molecular mechanisms leading to high-altitude pulmonary hypertension (HAPH) remains poorly understood. We previously analyzed the whole genome sequence of Kyrgyz highland population and identified eight genomic intervals having a potential role in HAPH. Tropomodulin 3 gene (TMOD3), which encodes a protein that binds and caps the pointed ends of actin filaments and inhibits cell migration, was one of the top candidates. Here we systematically sought additional evidence to validate the functional role of TMOD3. In-silico analysis reveals that some of the SNPs in HAPH associated genomic intervals were positioned in a regulatory region that could result in alternative splicing of TMOD3. In order to functionally validate the role of TMOD3 in HAPH, we exposed Tmod3^−/+ mice to 4 weeks of constant hypoxia, i.e. 10% O₂ and analyzed both functional (hemodynamic measurements) and structural (angiography) parameters related to HAPH. The hemodynamic measurements, such as right ventricular systolic pressure, a surrogate measure for pulmonary arterial systolic pressure, and right ventricular contractility (RV- ± dP/dt), increases with hypoxia did not separate between Tmod3^−/+ and control mice. Remarkably, there was a significant increase in the number of lung vascular branches and total length of pulmonary vascular branches (P < 0.001) in Tmod3^−/+ after 4 weeks of constant hypoxia as compared with controls. Notably, the Tmod3^−/+ endothelial cells migration was also significantly higher than that from the wild-type littermates. Our results indicate that, under chronic hypoxia, lower levels of Tmod3 play an important role in the maintenance or neo-vascularization of pulmonary arteries.

Introduction

High altitude (HA) environments exert a strong selection pressure leading to adaptation that involves genetic, epigenetic and physiological changes during the course of evolution. Although, adaptation primarily revolves around means to cope with the chronic hypoxic condition of HA, the physiological/pathophysiological response and the associated genes differ across species and between populations (1–7). For example, in the Andeans highlanders excessive erythrocytosis is common, and in the Tibetans right heart failure due to HA exposure is reported but without excessive erythrocytosis (8). An inevitable consequence of ascent to HA is the rise in pulmonary pressure or high-altitude pulmonary hypertension (HAPH) (9,10), as humans and animals who go to HA develop PH (11). Even though HAPH is the hallmark of most of the HA diseases, a consensus on its disease definition and pathogenesis is not yet attained (11,12). Since HAPH is prevalent in Kyrgyz highlanders (~14–20%) (13), we focused on this population and recently performed a whole genome sequence (wgs) analysis of Kyrgyz highlanders (2) in order to identify the DNA regions (and specific genes) that have been selected by evolution. Genome scans revealed strong selection signals in eight regions, encompassing multiple genes, with almost all of the single nucleotide polymorphisms (SNPs) in the noncoding or intergenic regions (2). This is rather consistent with most of the complex diseases or traits, where the phenotypes are driven by noncoding variants having regulatory functions, such as enhancers, alternate splicing and epigenetic changes (14).

We have proposed Tropomodulin-3 (TMOD3) as one of the candidate genes of HAPH in the Kyrgyz population (2). Since TMOD3 is highly expressed in lungs, we systematically sought additional evidence, some bioinformatic and others experimental (in vitro and in vivo), to validate the functional role of TMOD3. For the analysis, we used Genotype-Tissue Expression (GTEx) consortium data that comprehensively characterize genetic associations for gene expression and alternate splicing (15). In addition, we used Tmod3^−/+ mice to investigate experimentally the role of TMOD3 in hypoxia.

Results

High-altitude pulmonary hypertension associated SNPs overlay with alternative splicing/promoter sites of TMOD3 gene

Using a sliding window haplotype of 50 kbp, we identified eight genomic intervals under selection in Kyrgyz highlanders that constitutes several genes (n = 53) (2). In order to analyze the functional significance of the SNPs in the selected intervals, we examined their genomic positions with reference to the positions of various regulatory events, e.g. enhancers, splicing, methylation. SNPs overlaying such sites might have a role in gene activity, for example in relation to tissue specificity, and developmental activity. Upon systematic analysis of the genomic regions, we observed SNPs overlaying alternative promoter (AP) sites in one gene from each of the top two intervals (Supplementary Material, Figs S2 and S3). These are rs16943176 and rs12946397 overlaying AP for RAD51C from interval-1 (Supplementary Material, Fig. S2A and B) and rs2570233 and rs2570258 for TMOD3 from interval-2 (Fig. 1A, orange highlight, Supplementary Material, Fig. S3). Since APs are positively correlated with alternate splicing (AS) (16), we speculate that the SNPs overlaying AP sites affect the relative abundance of AS variants of the candidate genes, which in turn, could play key roles in the regulation of major biological processes, including response to environmental stress (17). In order to study the impact of the SNPs on endothelial TMOD3 level, we used two human umbilical vein endothelial cell (HUVEC) lines, each of A/A and A/G genotype, and exposed them to 1% O₂ for 24 h. Although, we could not obtain cells with G/G genotype, the RT-PCR results show a lower level of TMOD3 in the A/G cells when kept in 1% O₂ environment for 24 h (Supplementary Material, Fig. S11). Simultaneously, the ENST00000558455 level in A/G under both 21% O₂ and 1% O₂ were higher, indirectly indicating a higher functional TMOD3 level in the A/A genotype. These functional SNPs are in addition to our previously reported missense mutations in Heat Shock Transcription Factor 5 (HSF5, rs3803752 and rs117817367) (2). Interestingly, the tissue specific expression reveal both RAD51C and TMOD3 are similarly expressed (Supplementary Material, Fig. S4A).

Figure 1

SNPs associated with No-HAPH in the selected interval constitutes lung specific sQTL of TMOD3 gene. (A) SNPs (top line marks) near TMOD3 in the selected interval that differs significantly between HAPH versus No-HAPH, i.e. allele frequency differences are >50%. The ‘alternative event’ track shows the predicted alternative promoter. These are further highlighted/aligned with SNP position. The blocks represent exons; lines indicate introns, drawn such that no exons overlap. The bar plots represent the alternate allele frequency (%) distribution in the HAPH, No-HAPH and NES (method) of sQTL. (B) A representative violin plots measuring normalize intron-excision ratio of one of the SNPs (rs2570233) that is overlapping an alternative promoter (alternative splicing) sites (NES = −0.42). In the lungs, the SNPs are associated with the excision level of an intron (position chr15:51916133-51 934 250, identified as cluster clu_21347, Supplementary Material, Fig. S12), where the level of intron-defining reads (in the RNA-seq data of GTEx lung sample) flanking base position 51 916 133 to 51 934 250 on chromosome 15 is high in A/A, medium in A/G and low in G/G genotype (P-value < 0.0001). The clusters are identified using LeafCutter that employs split-reads to uncover alternative choices of intron excision by finding introns that share splice sites (Ref. 38). A detailed explanation on how LeafCutter identifies sQTL is provided in the Methods section. (C) The expression of TMOD3 significantly increases in HPAEC when cultured under 5% O₂ environment. Error bar represents SEM; ^*, P-value < 0.05.

HAPH associated SNPs are identified as splicing QTLs for TMOD3 in lung

In order to infer functional roles of these SNPs, we utilized GTEx Analysis Release V8 (dbGaP Accession phs000424.v8.p2) to identify splicing Quantitative Trait Loci (sQTL). GTEx RNA-seq data depicted 11 different AS isoforms for TMOD3, where the AS variant ENST00000308580 is the most abundant and functionally normal isoform, as it has the critical tropomyosin-binding sites, i.e. TMBS1 and TMBS2, and actin-binding sites, i.e. ABS1 and ABS2 coded by exons 3, 4, 6–10, 12 and 13 (Supplementary Material, Fig. S5A-C). Interestingly, other AS variants are un-translated or none of them retains all the critical exons (Supplementary Material, Fig. S5A-B), indicating these variants as functionally impaired. Additionally, the abundance of these AS variants also depend on different tissues suggesting a tissue specificity, e.g. ENST00000558666 is higher than ENST00000308580 in brain (Supplementary Material, Figs S5C and S6). Of note, GTEx dataset provides a comprehensive platform to analyze AS and identify sQTL as it has both short-read RNA-seq from different tissues/organs and their corresponding genotype information (11). Here we checked for significantly associated (P-value < 0.05) single-tissue sQTLs in the genes of selected intervals in all the tissues/organs and more specifically in lung. Out of the two genes having AP sites, i.e. RAD51C and TMOD3, no sQTL was detected for RAD51C. There were 835 SNPs (n = 410 for lung) significantly associated single-tissue sQTLs that cumulate to 8157 instances when combined from different tissues/organs for TMOD3 (Supplementary Material, Table S1). Interestingly, 170 of the 185 SNPs associated with HAPH were identified as lung-sQTLs for TMOD3, i.e. 140 of 410 SNPs (Supplementary Material, Tables S2 and S3). More importantly, all the SNPs of the interval-2 haplotype (n = 61, spanning 275 kbp) were among the significantly associated lung sQTL (Fig. 1A and B, Supplementary Material, Fig. S3). Nonetheless, it is important to note that the sQTLs are genetic variants associated with an intron-defining reads computed in a normalized space where magnitude may not have a direct biological meaning.

TMOD3 responds to hypoxia in lung

Our above results suggest that TMOD3 is a likely candidate gene of HAPH. To verify further that this candidate gene has an influence in the lung under hypoxia, we first examined its expression in GTEx organ samples and then investigated its response to a hypoxic environment, as at HA. The expression of TMOD3 is high in the lung as compared with other tissues in GTEx (Supplementary Material, Fig. S4B). In order to test its response to hypoxia, we treated HPAEC to 1% and 5% O₂ and quantified the temporal transcriptomic changes. Remarkably, the levels gradually but significantly increased (P-value < 0.05, Fig. 1C and Supplementary Material, Fig. S7). The increase is significant as early as 2 h and is maintained continuously to 48 h after exposure (Fig. 1C and Supplementary Material, Fig. S7).

Hemodynamic changes in Tmod3^+/+ and Tmod3^−/+ after chronic hypoxia

The relative transcript levels of Tmod3 in Tmod3^−/+ mice were half of the levels in Tmod3^+/+ mice both at the tissue level, i.e. Lung and Spleen as well as in the endothelial cells, confirming the functional impact of the heterozygous genotype (Supplementary Material, Fig. 1B). In order to investigate the role of TMOD3 in hypoxia-induced pulmonary hypertension, we measured hemodynamics parameters in Tmod3^+/+ and Tmod3^−/+ mice exposed to 10% O₂ for 4 weeks. As anticipated, a significant rise in RVSP, mean pulmonary artery pressure (mPAP) and RV- ± dP/dt in the Tmod3^+/+ mice occurred after such chronic hypoxia exposure (Fig. 2). The mean RVSP in Tmod3^+/+ mice increased from 22.1 ± 0.9 mmHg to 32.2 ± 2.3 (Fig. 2C). Similarly, the mPAP and RV-dp/dt increased from 15.3 ± 0.5 and 1882 ± 124 to 21.3 ± 1.4 and 2669 ± 689 mmHg/s, respectively (Fig. 2D and E). A similar increase for RVSP, mPAP and RV- ± dP/dt was also detected in the heterozygote mice Tmod3^−/+ (Fig. 2) with no significant difference the wild-type Tmod3^+/+ mice (Fig. 2D and E). The comparisons of fulton index between Tmod3^−/+ and Tmod3^+/+ mice measured after 4 weeks of hypoxia treatment remain non-significant (P-value = 0.44, Supplementary Material, Fig. S8). Interestingly, heart rate (HR) remained in a similar range under both environmental conditions, i.e. 21% and 10% O₂ exposure and for both Tmod3^+/+ and Tmod3^−/+ mice (Fig. 2F).

Figure 2

Hemodynamic changes in Tmod3^+/+ and Tmod3^−/+ after treatment with 4 weeks (chronic) hypoxia are comparable. (A) Representative graphs of RVP and (B) RV ± dP/dt. (C) RVSP; (D) mPAP and (E) RV ± dp/dt_max measurements depict significant increase under hypoxic, i.e. 10% O₂, condition in both Tmod3^+/+ and Tmod3^−/+ mice. (F) HR remain comparable in 21% and 10% O₂ environments. No differences were detected between the Tmod3^+/+ and Tmod3^−/+ mice when compared under respective O₂ environments. Normoxia, Tmod3^+/+ (n = 4) and Tmod3^−/+ mice (n = 4); hypoxia, Tmod3^+/+ (n = 6) and Tmod3^−/+ mice (n = 9). Error bar represent SEM.

Figure 3

Improved vascularization in the lungs of Tmod3^−/+ mice after 4 weeks of chronic hypoxia treatment. (A) Representative angiography images of the lungs of Tmod3^+/+ and Tmod3^−/+ mice under normoxia and after 4 weeks of 10% O₂ treatment. Image of whole lung at 8x and part of the lung at 30x are simultaneously depicted. (B) The total length of branches (mm/mm²) increase significantly in the Tmod3^−/+ mice at 4 weeks of 10% O₂ exposure. (C) Number of branches per mm² also depicts a significant increase in the Tmod3^−/+ mice under hypoxic environment. (D) Chronic hypoxia exposure also indicates an increasing trend in the number of junctions per mm². Normoxia, Tmod3^+/+ (n = 3) and Tmod3^−/+ mice (n = 4); hypoxia, Tmod3^+/+ (n = 3) and Tmod3^−/+ mice (n = 5). ^*, P-value < 0.05 (Mann–Whitney test).

Better vascularization in the lungs of Tmod3^−/+ mice after chronic hypoxia

Under normoxic condition, the vascular density in the Tmod3^−/+ mice is similar to that in control mice. One of the pathophysiological consequences of chronic hypoxic exposure, especially in rodent models, is vascular pruning, i.e. loss of the peripheral pulmonary vasculature shown in angiogram and impaired vascular regeneration. Images in the control mice (Tmod3^+/+) revealed a significant decrease in the number of lung vascular branches under 10% hypoxia (Fig. 3A-D). For instance, the number of pulmonary vascular branches per unit area (mm²) decreased from 620.4 ± 67.5 to 201.6 ± 89.5 (Fig. 3C; P-value < 0.05) and the total length of pulmonary vascular branches decreased from 22.7 ± 1.0 to 15 ± 1.3 mm/mm² (Fig. 3B; P-value < 0.05). The number of vascular junctions also decreased, although not significantly, from 300.0 ± 80.4 to 155.9 ± 114.6/mm² (Fig. 3D; P-value = 0.20). However, remarkably, under 10% O₂, the pulmonary vascular density was significantly higher in the Tmod3^−/+ mice when compared with the Tmod3^+/+ mice, as indicated by the number of lung vascular branches, the number of lung vascular junctions and the total length of lung vascular branches (Fig. 3). For example, the number of branches under 10% O₂ was 470.4 ± 46.9 per mm² which is close to the density under normoxia at 570.8 ± 128.8 per mm² (Fig. 3C). Similarly, the number of junctions and the total lengths of branches also remained similar, i.e. 335.3 ± 40.0 per mm² and 21.8 ± 2.3 mm/mm² in 21% O₂ and 293.4 ± 42.0 per mm² and 21.5 ± 1.1 mm/mm² in 10% O₂, respectively (Fig. 3B and D).

To examine any difference in the hypoxia sensitivity of the pulmonary vessels of Tmod3^−/+ and control mice, and to rule out the possibility that the increased vascularization or lack of decrease vascularization in the 4-week hypoxia treated Tmod3^−/+ mice was not due to any lack of vessel ‘pruning’ in the Tmod3^−/+ mice we examined the basal level and rise in PAP in an ex-vivo analysis of lungs obtained from normoxia and 4-week hypoxia-treated Tmod3^+/+ and Tmod3^−/+ mice (methods). The basal PAP in the isolated perfused/ventilated lung preparation increased significantly in the lungs obtained from the mice kept under 10% O₂ environment irrespective of their genotype (Supplementary Material, Fig. S9A and B). Although the hypoxia-induced rise in PAP is similar in the lungs of mice kept under both 21% and 10% O₂ environment (Supplementary Material, Fig. S9A and C), the 40 mM K⁺-induced increase in PAP, depicting membrane depolarization-induced pulmonary vasoconstriction, is lower in the mice kept in 10% O₂ (Supplementary Material, Fig. S9A and D).

Faster Tmod3^−/+ endothelial cell migration

In order to further understand the function of Tmod3 and its role in HAPH, we studied in vitro the differences between wild type Tmod3^+/+ and Tmod3^−/+ in endothelial cell migration. We isolated endothelial cells from Tmod3^+/+ and Tmod3^−/+ mice and performed scrape/injury assay under hypoxic conditions (methods). Remarkably, the results indicated a significantly faster increase in the Tmod3^−/+ endothelial cell migration than in the wild type (Fig. 4). Cell migration was measured at 6, 12, 18 and 24 h after the initial scrape/injury. The percentage of the area covered/healed by the endothelial cells reflected the migration of these cells and our analysis indicated an increase in Tmod3^−/+ cell migration at 6 h, i.e. 26.1% versus 33.6% with statistically significant differences after 12 h, i.e. 55.6% versus 72.4% at 12 h (Fig. 4B). At 18 h, the scrape/injury area of the Tmod3^−/+ cells was almost completely covered (>95%) compared with 82% in the Tmod3^+/+ cell (Fig. 4B).

Figure 4

Increase cell migration of Tmod3^−/+ endothelial cells. (A) Representative image of endothelial cell migration for cells isolated from the lungs of Tmod3^+/+ and Tmod3^−/+ mice. (B) Area covered by endothelial cells at 6, 12 and 18 h after scrape/injury showed that the cells isolated from Tmod3^−/+ mice have significantly increased restorative ability (migration) compared with the endothelial cells isolated from control mice, i.e. Tmod3^+/+. P-values, ^* = 0.0178 and ^** = 0.0078. Error bar represents SEM.

Discussion

The present study examined the impact of the gene Tmod3 on structural, functional and hemodynamic changes in lung under chronic hypoxia treatment. We have demonstrated that the lung of the heterozygote Tmod3^−/+ mice developed a significant increase in vascular density when exposed to 10% O₂ for 4 weeks (Fig. 3). Since the loss of precapillary vessels in the lungs (pruning) and impaired vascular regeneration have been well established in hypoxia-induced PH (18,19), a higher lung vascularization in the hypoxia treated Tmod3^−/+ mice (Fig. 3) is indeed remarkable and has direct implications on HAPH. Loss of precapillary vessels in chronic hypoxia exposure is due to a significant death and abnormal proliferation and migration of the pulmonary endothelial cell (20,21) which, along with fibroblast and smooth muscle cells proliferation, leads to abnormal vascular remodeling (22). However, it is important to understand whether the higher vascular density in the Tmod3^−/+ mice is due to neovascularization or to diminished hypoxia sensitivity/response. In humans, the tropomodulin (Tmod) family includes four homologs (TMOD1–4), all having the role of pointed-end capping of actin filament regulating actin filament architecture in diverse cell types (23,24). However, different Tmod homologs exhibit some level of tissue specificity. For example, TMOD1 is expressed in striated muscle, red blood cells, lens fiber cells and neurons, TMOD2 is neuron specific and TMOD4 is predominantly expressed in the skeletal muscle. TMOD3 is although ubiquitously expressed but its expression is relatively higher in the lungs (25,26). Given the diverse cell types the TMOD3 is expressed, recent studies have reported its involvement in determination of cell shape, endothelial cell migration, muscle contraction etc. (26). It is important to note, in this regard, that TMOD3 is upregulated under hypoxic condition (Fig. 1C and Supplementary Material, Fig. S7) in the endothelial cells, presumably due to TMOD3 being an Akt2-interacting partner (27). This increased level may indeed significantly hinder endothelial cell migration during chronic hypoxia, resulting into abnormal proliferation and migration (24,28). Suppression of TMOD3, as in Tmod3^−/+, increases endothelial cell motility, thus leading to vascular regeneration (Fig. 3). This is actually quite similar to what occurs in proliferative retinopathy, a pathological retinal angiogenesis (28,29). In addition, an alteration of TMOD3 in the lungs may have more implication since TMOD3 transcript level is higher in the lungs (Supplementary Material, Fig. S4B) and it is the only Tmod homolog expressed in the endothelial cells, unlike in other tissues (28,30). Therefore, in the motile endothelial cells, especially during hypoxia-induced vascular remodeling in the lungs, a lower TMOD3 would dysregulate the dynamics of actin filament lengths during lamellipodia extension, thus promoting cell migration (31).

The unexpected finding that depicts no change in the Tmod3^−/+ hemodynamics after hypoxia in the presence of greater vascularization than in control lungs is intriguing. The expected inverse relationship between lung vascularization and vascular resistance did not occur in Tmod3^−/+. The ambiguity may arise due to a higher than expected resistance or flow. Since the RV- ± dP/dt_max (mmHg/s) does not differ between Tmod3^−/+ and Tmod3^+/+, a possible explanation is that Tmod genes are essential for the maintenance of contractile actomyosin bundles in muscle cells (31,32) and hence may play a role in smooth muscle contraction. Indeed, in addition to its primary role of capping pointed end (−) of the actin filament and inhibiting polymerization and depolymerization of actin monomers (24,33,34), Tmod3 also regulates sarcoplasmic reticulum (SR)-associated γ-actin architecture, mechanically stabilizing the sarcoplasmic reticulum via a cytoskeletal linkage to sAnk1.5. (23). Since SR serves as a Ca²⁺ reservoir, critical for vascular smooth muscle contraction, a downregulation of Tmod3 (Tmod3^−/+) will not only disrupt the dynamic regulation of thin filament’s length (especially F-actin stability) (35) but can also cause SR morphological defects leading to impaired Ca²⁺ release (23) and vascular smooth muscle contraction, a potential explanation for the lack of any change in the PAP in the hypoxia-treated Tmod3^−/+ mice. Furthermore, a recent study shows that Tmod genes control the balance between protrusive and contractile structures by stabilizing actin–tropomyosin filaments (30). This study specifically shows that Tmod3 along with Tmod1 are essential genes for the maintenance of contractile actomyosin bundles and that Tmod-dependent capping of actin–tropomyosin filament as critical for the regulation of actin homeostasis in non-muscle cells (30).

Since TMOD3 is one of the genes in the top selected interval of the genomes of Kyrgyz highlanders, the current study in Tmod3^−/+ mice provides functional evidence demonstrating that a lower level of TMOD3 may have a critical role in potentially protecting no-HAPH Kyrgyz highlanders from developing severe PH. The hypothesis of a lower TMOD3 in the no-HAPH individuals proposed here is supported by a higher alternative allele frequency for the SNP variants overlaying AP sites of TMOD3 that are also identified as sQTL in the lung (Fig. 1). Since the coordination between tissue-specific AP is mediated at least in part by tissue-specific regulatory factors with conserved sequence motif in APs (36), it is likely that, in no-HAPH, a higher frequency of alternate alleles may result in an unusual AP sequence motif for cell-specific expression of normal TMOD3, and therefore a lower TMOD3 lung-specific transcript level. Additionally, abundance of any other splice variant will result in a dysfunctional TMOD3 isoform, as would be expected because of the absence of TMBSs and ABSs regions (Supplementary Material, Fig. S5). Finally, with all considered, we do agree that the core issue of HAPH pathogenesis, i.e. hypoxia-induced rise in PAP, still remains elusive as the PAP in Tmod3^−/+ mice remains high. We believe that the increased vascularization, as seen in the Tmod3^−/+ mice in hypoxia, is the result of only a modest proportion of genome selection to HA among Kyrgyz highlanders. Hence, we believe that there are additional genes embedded in the eight genomic intervals that we identified (2) that have still a major role to play in order to understand the full pathophysiology of no-HAPH adaptation.

There are certain limitations to our current observation. First, our HAPH analysis is limited to human wgs analysis and the RNA/proteins could not be obtained for validation because of multiple reasons, including population access of samples. However, based on the data from the GTEx Consortium atlas of regulatory genetic variations which is based on genotype-RNASeq data and depicts its effects across human tissues (15), the SNPs we identified as associated with HAPH is also associated in a dose-dependent manner with the splice variants (sQTLs) of TMOD3. This leads us to speculate that a SNP frequency difference of > 50% between cases and controls is partly contributing to the HAPH phenotype in the Kyrgyz population. In conclusion, the results presented here, together with our previous findings from wgs analysis of Kyrgyz population (2) reveals that a lower level of primary TMOD3 transcript likely associates with No-HAPH in the Kyrgyz highlanders. As we understand the diverse functional role of different Tmods, we can potentially delay structural remodeling or restore existing precapillary vascular bed by utilizing prospective Tmod3 blockers that would substantially delay PH progression. More importantly, our result reveals a fine example of the critical link between molecular and evolutionary biology.

Materials and Methods

Analysis of the individual SNPs and their plausible regulatory function

Previously, we identified selection signals in eight intervals, encompassing multiple genes, associated with the phenotype of the adapted subjects in the Kyrgyz highlanders (No-HAPH) (2). In the current study, our endeavor was to evaluate TMOD3 as a specific candidate gene for having a role in HAPH. For this, we systematically evaluated first the individual SNPs of the selected intervals. This includes looking into the genomic positions of the SNPs to see if they align with any regulatory elements. We referred to the annotations in ‘UCSC Alt Events Track Settings,’ ‘SIB Alt-Splicing Track Settings,’ on the UCSC genome browser version Human GRCh37/hg19. Simultaneously, we also checked if the SNPs are identified as eQTL or sQTL, specifically in the lungs, using GTEx portal (15). We begin our analysis from the most significant interval, i.e. interval-1. On failing to detect any evidence that can direct us to a specific candidate gene, we moved on to the next most significant interval in line.

sQTL analysis

To identify sQTL that significantly associate with No-HAPH, we utilize GTEx analysis release V8 (dbGaP Accession phs000424.v8.p2) that has both short-read RNA-seq and their corresponding genotype data (15). A detailed method for sQTL analysis is available on https://gtexportal.org/home/documentationPage. The genotype data used for sQTL analyses in GTEx portal were based on wgs from 838 donors, which all had RNA-seq data available in V8. For sQTL analysis, GTEx portal included the variants with minor allele frequency ≥ 1% across all 838 samples. The cis-sQTL mapping was performed using FastQTL (37). As indicated, splicing is quantified using the intron excision phenotypes computed by LeafCutter (38). Briefly, LeafCutter utilizes split-read to construct a graph where nodes are introns and edges represent shared splice junctions between two introns. The connected components define intron clusters (as depicted in Fig. 1a in Li et al. [38]), namely that the LeafCutter obtains read proportions for all introns within alternatively excised intron clusters. The values are then standardized across individuals for each intron and quantile normalized across introns and used as the phenotype matrix. The linear regression is used (as implemented in fastqtl) to test for associations between variants (MAF ≥ 0:05) within 100 kb of intron clusters and the rows of the phenotype matrix that correspond to the introns within each cluster (Li et al. [38]).

In order to test whether a SNP is single-tissue sQTLs, nominal P-values are generated for each SNP by testing association between genotype and gene expression. The mapping window was defined as 1 megabase up- and down-stream of the transcription start site. All splicing phenotypes mapping to a gene were mapped jointly, using grouped permutations (-grp option in FastQTL). Beta distribution-adjusted empirical P-values from FastQTL were used to calculate q-values (39), and a false discovery rate (FDR) threshold of ≤ 0.05 was applied to identify genes with a significant sQTL. The normalized effect size (NES) of the sQTLs is the slope of the linear regression, and is computed as the effect of the alternative allele relative to the reference allele in the human genome reference GRCh38/hg38 (i.e. the sQTL effect the alternative allele).

Mice

Tmod3^−/+ heterozygous mice were generated by Dr Velia M. Fowler, The Scripps Research Institute, La Jolla, CA (19). Since the Tmod3^−/− KO mice are embryonic lethal (40), we used Tmod3^−/+ as breeders. DNA was isolated from the tail prep using REDExtract-N-Amp™ Tissue PCR Kit (Sigma-Aldrich). Tmod3^+/+ and Tmod3^−/+ genotypes were confirmed using primers previously reported by Sui et al. (40). Briefly, a PCR product size of 729 bp band depicts wild type Tmod3^+/+ and an additional band size of 313 bp indicates Tmod3^−/+ mice (Supplementary Material, Fig. S1A). Both Tmod3^+/+ and Tmod3^−/+ were phenotypically normal under ambient conditions.

Hypoxia exposure and hemodynamic measurements

In order to rule out any hormonal-related influence on pulmonary hypertension during the estrous cycle (41), all animal experiments were conducted in male mice (8–10 weeks old, ~25 g body weight). Animal care and handling was performed according to the NIH guidelines. The protocols were approved by the Institutional Animal Care and Use Committee (IACUC) of the University of California, San Diego (UCSD).

For the chronic hypoxia treatment, mice were exposed to normobaric hypoxia where the O₂ was constantly maintained at 10% in a well-ventilated chamber for four weeks. In order to measure right ventricle pressure (RVP) and RV contractility, we catheterized right heart using a pressure catheter (Millar Instruments, PVR1030, 1 F, 4 E, 3 mm, 4.5 cm, Colorado, USA) introduced via the external right jugular vein. During the procedure mice were continuously anesthetized under 1.5% isoflurane. Data were recorded and analyzed using Lab Chart Pro1.0 software (ad Instruments). To evaluate RV hypertrophy, we isolated the whole heart, and then we separate the RV from the left ventricle (LV) and septum (S), and measured the weight of the RV, LV and S. Fulton Index, the ratio of the weight of RV to the weight of LV and S [RV/(LV + S)) was calculated to determine RV hypertrophy.

Measurement of hypoxic pulmonary vasoconstriction in ex vivo model, isolated perfused/ventilated mouse lung

Tmod3^+/+ and Tmod3^−/+ mice were anesthetized by pentobarbital sodium (120 mg/kg) via intraperitoneal injection. A detailed method is described by Jain et al. (42,43). Briefly, after tracheostomy, isolated lungs were immediately ventilated with normoxic gas mixture of 21% O₂/5% CO₂ using a rodent ventilator (Minivent type 845, Harvard Apparatus, USA). Respiratory rate and tidal volume were maintained at 80 breaths/min and approximately 250 μL. A 2 cmH₂O of positive end expiratory pressure was maintained. A pressure transducer (MPX type 399/2, Hugo Sachs Elektronik-Harvard Apparatus, Germany) connected to a tracheal catheter was used to measure the end inspiratory pressure. Anesthetized mice were placed in an isolated lung open perfusion system chamber (IL-1 Type 839, Harvard Apparatus, USA) with a heated water jacket at 37°C. After tracheal intubation, the chest was opened by median sternotomy and thymus and adipose tissue were carefully excised. Heparin (20 IU) was immediately injected into the right ventricle to prevent blood from coagulation.

To continuously measure pulmonary arterial pressure (PAP), a catheter connected to pressure sensor (P75 Type 379, Hugo Sachs Elektronik-Harvard Apparatus, Germany) was inserted into the main PA via the right ventricle. Another catheter was inserted into the left atrium via a small incision of the LV to allow perfusate to drain to reservoir. The pulmonary flow rate was set and maintained at 1 mL/min by a peristaltic pump (ISM 834, ISOMATEC, USA). The Powerlab data acquisition system (ad Instruments, CO, USA) was used to store and analyze the imaging data.

The lung vasculature was consistently superfused with Physiologic Salt Solution (PSS) composed of 120 mM NaCl, 4.3 mM KCl, 1.8 mM CaCl₂, 1.2 mM MgCl₂, 19 mM NaHCO₃, 1.1 mM KH₂PO₄, 10 mM glucose and 20% fetal bovine serum (pH 7.4) via a pump, whereas the lung airways and alveoli were ventilated with normoxic or hypoxic gas. Raising extracellular [K⁺] from 4.7 to 40 mM in PSS causes membrane depolarization in PASMC and pulmonary vasoconstriction due to a shift of the equilibrium potential for K⁺ from −85 to −31 mV. Before experimentation, the isolated lungs were first superfused with the 40 mM K⁺-containing PSS (40 K), at least three times, to stabilize the basal PAP and the amplitudes of 40 K-induced increases in PAP. When the basal PAP was stabilized, the lungs were repetitively challenged by ventilation of hypoxic gas (1% O₂ in N₂, for 4 min) to induce an increase in PAP due to alveolar hypoxia-induced pulmonary vasoconstriction (HPV). In the interval of hypoxic challenges, the lungs were ventilated with normoxic gas (21% O₂ in N₂).

Lung angiography

A detailed method on lung angiography is available in Jain et al. (42). Briefly, the anesthetized mice (using intraperitoneal injection of pentobarbital sodium, 120 mg/kg) were immediately injected with heparin (20 IU) into the heart to prevent blood from clotting. The lung was perfused with phosphate buffer saline via a polyethylene (PE-20) tube cannulated into the PA via the right ventricle using an automated pump (NE-300, Pump Systems, for 3 min at a speed of 0.05 mL/min). Subsequently, 0.08 mL of microfil polymer (yellow) (FlowTech Inc., Carver, MA) was perfused into the PA at a speed of 0.05 mL/min. The microfil polymer-filled lungs were kept at 4°C overnight. The next day, the lungs were dehydrated using different concentrations of ethanol: once in 50%, 70%, 80% and 95% ethanol, and twice in 99.9% ethanol. After dehydration, the lungs were placed in methyl salicylate (Sigma Aldrich, USA) at room temperature on a shaker for overnight in order to show only the lung vasculature. Lungs were then photographed or imaged with a digital camera (MU1000, FMA050, Amscope, CA). The peripheral lung vascular image, covering the peripheral area of the lung, 1 mm width from the edge was selected with Photoshop CS software, and the branches on the images in Photoshop were outlined manually and later converted to binary images with NIH Image J 1.8v software for quantitative analysis. The total length of branches, the number of branches and the number of junctions on the skeletonized images were obtained by Image J software and were normalized by the area selected within the peripheral regions of the lung.

Cell culture

Primary Pulmonary Artery Endothelial Cells (HPAEC) was purchased from ATCC (PCS-100-022™). For maintaining normal growth, we followed protocol as indicated by ATCC. The Endothelial Cell Growth Kit-BBE (ATCC® PCS-100-040) was added as indicated. We passaged the cells when cultures reached ~80% confluence. For the hypoxia experiments, equal numbers of cells were plated in five 60 mm cell culture dishes and maintained in regular incubator of room air, 5% CO₂ and 37°C. On day 3 of cell expansion, i.e. ~ 70% confluence, four plates were transferred to an incubator with 1% O₂, 5% CO₂ and 37°C. Cells from each dish were used to isolate RNA at 0, 2, 8, 24 and 48 h of hypoxia exposure. The experiments were conducted in triplicate and each one is from different passage.

Endothelial cell isolation and culture

For in-vitro experiment in mice, we used a modified protocol for isolating mouse lung endothelial cells (44). Briefly, the ICAM2 and CD31 antibody-coated dynabeads were utilized for endothelial cell enrichment. Mouse lung tissue dissociation and digestion we used collagenase I at 37°C for 45 min. This was followed by first cell sorting using CD31-coated beads and second sorting using ICAM2-coated beads. Culture media was changed every 48 h. Cells are ready for experiments when reach confluent monolayer. The cells were immediately used for experiments without further passage after second sorting.

Genotyping

DNA was isolated from eight cell lines (four iPSCs lines, two primary cell lines (one each of human primary artery endothelial cells and pulmonary artery smooth muscle cells) and two human umbilical vein endothelial cells (HUVEC)). The idea for including iPSCs was to reprogram these cells into endothelial cells if they were of the required genotype. For genotyping rs2570233 PCR and sequencing was done using specific primers TMOD3_seq_F: GCTTCGGGGAAGTCCATAGT and TMOD3_seq_R: CTCCCAAAGAGTTGCCAGAC and the sequencing results were viewed on Sequence Scanner Software V2.

Real-time qRT-PCR

RNA was isolated from tissue samples using RNeasy Mini Kit (Qiagen, USA). We used SYBR® Green Master Mix for RT-PCR which is a pre-formulated, optimized, universal 2X master mix for real-time PCR workflows. RT-PCR was performed on CFX96 Real-Time PCR System. The specific primers used for RT-PCR includes TMOD3_F: TCTGGAGAAAGAAGCATTGGA and.

TMOD3_R: CAGGTTTCTGTTTGGGGATAA for humans and Tmod3_F2: CAAGCATTGGAGCACAAAGA and Tmod3_R2: ACATTGGGAAAACGCTCTTG for mice (45). Additionally, the specific primers used for RT-PCR of splice variant ENST00000558455 were TMOD3_558455_F: GTTTCCGTGGGTCTGTGAAG and TMOD3_558455_R: GCCCTGGATCAATCTCAGTC. The relative transcript levels of Tmod3 in Tmod3^−/+ mice were compared with its level in Tmod3^+/+ mice in the Lung and Spleen as well as in the endothelial cells isolated from lungs.

Endothelial cell migration assay

Endothelial cells from Tmod3^+/+ and Tmod3^−/+ were cultured in six well plates. Cell monolayers were wounded by scraping with a pipette tip and incubated at 1% O₂, 37°C for 24 h. Images of the cells were taken immediately after scraping and after every 6 h to measure cell migration into wounded areas. The area covered by the migrating cells was quantified and normalized to the original wound area (n = 12 per group).

Statistical analysis

The composite data are shown as means ± standard errors (SEM). Paired or unpaired Student’s t-test and one-way analysis of variance (ANOVA) with Bonferroni multiple comparison test were used for statistical analysis. A P-value < 0.05 was considered as statistically significant.

Acknowledgements

We thank Prof. Velia M. Fowler, The Scripps Research Institute, La Jolla, CA for providing the Tmod3^-/+ mice. We thank Dr Wei Wu and Dr Huiwen Zhao for providing the HUVEC and iPSC lines. We thank Travis Smith, Juan Wang and Orit Poulsen for technical support.

Conflict of Interest Statement. No conflicts of interest to declare.

Funding

National Institutes of Health (R01HL127403 to G.G.H., and R01GM114362 to V.B.).