Rapid comparison of surface protein isoform membrane topologies by surfaltr

surfaltr Installation

surfaltr is freely available for download via Bioconductor. surfaltr has some dependencies which should all be installed/imported automatically when the code shown below is run to install surfaltr:

if (!require("BiocManager", quietly = TRUE))
    install.packages("BiocManager")

BiocManager::install("surfaltr")

To load surfaltr, the following code can be used:

library(surfaltr)

The dependencies of the package that are automatically imported include: dplyr (>= 1.0.6), biomaRt(>= 2.46.0), protr (>= 1.6-2), seqinr (>= 4.2-5), ggplot2 (>= 3.3.2), utils (>= 2.10.1), stringr (>= 1.4.0), Biostrings(>= 2.58.0), readr (>= 1.4.0), httr (>= 1.4.2), testthat (>= 3.0.0), xml2 (>= 1.3.2), msa(>= 1.22.0), methods (>= 4.0.3)

The dependency of the package which is NOT automatically imported: TMHMM standalone software

surfaltr predicts TM domains either using TMHMM standalone software OR using the Phobius API. If the user intends to use the TMHMM fuctionality of surfaltr then they must install the TMHMM standalone software (details below).

When to Use Phobius or TMHMM

In terms of prediction specificity, Phobius’s ability to detect signal peptides contributes to a greater amount of detail than TMHMM. Additionally, as all sequences are run using an API as opposed to local software, the run time for Phobius is also slightly lower. However, if a protein has no TM domains and no signal peptide, Phobius automatically flags the protein as internal, which is not always true and should be verified experimentally before drawing any conclusions. If sequences are considered private information and should not be posted externally, we recommend using TMHMM to ensure everything is run locally.

Phobius Installation

As run_phobius() relies on the Phobius API, a copy of the software does not need to be downloaded on the user’s device and will instead automatically be called by the surfaltr functions.

TMHMM Standalone Software Installation

In order to be able to use TMHMM⁴ within surfaltr to predict membrane topology, it is important to first ensure that you have TMHMM 2.0 standalone software installed on your computer. To do this,
1. Navigate to https://services.healthtech.dtu.dk/service.php?TMHMM-2.0   2. follow directions for installation of Linux version of standalone software (tmhmm-2.0c.Linux.tar.gz).
3. Untar the tmhmm-2.0c.Linux.tar.gz file.
4. Extract it to a place where you store software.
5. Insert the correct path for perl 5.x in the first line of the scripts bin/tmhmm and bin/tmhmmformat.pl (if not /usr/local/bin/perl or #!/usr/bin/env perl). Use which perl and 6. perl -v in the terminal to help find the correct path.
6. Make sure you have an executable version of decodeanhmm in the bin directory.
7. Make sure to record the path to which TMHMM 2.0 is installed, as this information will be needed in order to run TMHMM in future functions.
TMHMM license and thus the TMHMM 2.0 standalone software is freely available to all academic institutions upon request from the official website linked above.

Troubleshooting TMHMM Installation

If the TMHMM standlone software is NOT installed in the correct path or folder you will receive error messages when running surfaltr functions. In this scenario please check the directory in which TMHMM2.0c was downloaded and make sure that the directory structure is as such: path/to/folder/TMHMM2.0c/bin/decodeanhmm.linux

TMHMM may get downloaded with the following directory names and structure: path/to/folder/TMHMM2.0c/tmhmm-2.0c/decodeanhmm.Linux_x86_64 Please fix the directory structure to be the one shown above to use TMHMM within surfaltr. OR run the following code:

check_tmhmm_install(path/to/folder/TMHMM2.0c)

This function also prints a helpful method providing tips on how to fix the installation if TMHMM is not found at the folder path specified. If TMHMM installation is problematic for some reason we recommend only using the Phobius method within surfaltr to predict TM domains.

Why filter genes known to encode surface proteins?

While surfaltr can be used to plot topology of several proteins with the limit being set by the bulk processing power of TMHMM standalone software and Phobius web server, we recommend limiting input to isoforms of genes known to encode surface proteins. This can be done by using existing protein repositories and annotation systems such as Uniprot and Surfacegenie⁶ to filter out isoforms that map to genes known not to encode surface proteins. This pre-processing step will ensure optimum performance of the package.

Gene name and transcript ID’s input (Input type 1)

The example input table containing data from Uhlén et al.⁷ is shown below with column names “gene_name” and “transcript”. This matches the format required for surfaltr to be able to process ID-based input data. This exmaple dataset is provided along with the package for the user and can be loaded

Gene name and amino acid sequence input (Input type 2)

The example input table containing data from Ray et al.⁸ is shown below with column names “external_gene_name”, “transcript_id”, and “protein_sequence”. This matches the format required for surfaltr to be able to process amino acid sequence-based input data.

From known Ensembl Transcript Models (Input type 1)

When a .csv file containing known Ensemble transcript IDs and gene names is provided, the graph_from_ids() function can be utilized to obtain the completed plot of paired isoforms. Additionally, the function can include optional inputs summarized in the table below:

An example of the code required to run this function applied to the Uhlén et al.⁷. data is shown below:

graph_from_ids(data_file = system.file("extdata", "hpa_genes.csv",package="surfaltr"), organism = "human", 
               rank = "combo", n_prts = 20, mode = "phobius", size_txt = 2, 
               space_left = -400, tmhmm_folder_name = NULL)

From Amino Acid Sequence Input (Input type 2)

When a .csv file containing gene names and amino acid sequences is provided, the graph_from_aas() function can be utilized to obtain the completed plot of paired isoforms. Additionally, the function can also include the optional inputs summarized in the table above.

An example of the code required to run this function applied to the Ray et al.⁸ data is shown below:

graph_from_aas(data_file = system.file("extdata", "CRB1.csv", package="surfaltr"), organism = "mouse", rank = "TM", 
               n_prts = 10, mode = "phobius", size_txt = 2, space_left = -400, 
               tmhmm_folder_name = NULL)

Output 1: Paired Isoforms

Below is a table depicting the table output from the get_pairs function applied to the Ray et al.⁸ dataset:

Output 2: FASTA File

Below is a snippet of the FASTA file produced from this function when applied to the Ray et al.⁸ dataset:

Using TMHMM

In order predict TM topology using TMHMM 2.0⁴, surfaltr includes a function called get_tmhmm(). This function uses the AA.fasta file generated as the output of the get_pairs() function and the full path to the folder containing the TMHMM 2.0 software as input (fasta_file_name). The path to TMHMM 2.0 should be provided as a string and end in “TMHMM2.0c”. For the fasta file input parameter, the input should be provided as a string that just includes the file name of the fasta file to be used without a path. The function will return a data frame (output as path to/working directory/mem_topo_tmhmm.csv) with the TMHMM-predicted location for each amino acid in all the protein sequences. The topology is expressed as i- inside the plasma membrane, M- TM, o- outside the plasma membrane and O – unknown/assumed to be outside.

To run get_tmhmm alone, the following code can be used:

mem_topo <- get_tmhmm(fasta_file_name = "AA.fasta", tmhmm_folder_name = "~/path/to/TMHMM2.0c")

It’s important to note that in order for this function to run, there needs to be a FASTA file containing the amino acid sequences for each transcript within the “output” folder in the working directory. If the get_pairs function is run before the get_tmhmm() function, this folder and file will be created automatically. The code below can be used to automate the fasta file generation and membrane topology prediction process:

AA_lst <- get_pairs(data_file = system.file("extdata", "CRB1.csv", package="surfaltr"), if_aa = TRUE, organism = "mouse")
mem_topo <- get_tmhmm("AA.fasta", tmhmm_folder_name = "~/path/to/TMHMM2.0c")

When to Use TMHMM

If the transcripts being run are not public information or should be kept private, the TMHMM standalone software along with the tmhmm R package provides an avenue to evaluate TM topology without using an external server, making it the ideal TM prediction algorithm in this instance. However, TMHMM is not able to predict the presence of signal peptides, possibly resulting in inaccuracies in TM prediction at the beginning of the sequence. Additionally, if a protein has no TM domains, TMHMM automatically flags the protein as extracellular, which is not always true and should be either verified experimentally or verified by comparison to the prediction by Phobius.

Using Phobius

In order predict TM topology using Phobius, surfaltr includes a function called run_phobius(). This function takes in the AA.fasta file generated as the output of the get_pairs() function (fasta_file_name). Additionally, this function is reliant on the data frame returned by the get_pairs() function as an input parameter (AA_seq). The function will return a data frame (output as path to/working directory/mem_topo_phobius.csv) with the Phobius-predicted location for each amino acid in all the protein sequences. The topology is expressed as i- inside the plasma membrane, M- TM, o- outside the plasma membrane and s – signal peptide.

To run run_phobius alone, the following code can be used:

topo <- run_phobius(AA_seq = AA_lst, fasta_file_name = "~/<your-workin-directory>/output/AA.fasta")

To automatically generate the FASTA file used by the function, the following code can be used:

AA_lst <- get_pairs(data_file = "~/CRB1.csv", if_aa = TRUE, organism = "mouse")
mem_topo <- run_phobius(AA_seq = AA_lst, fasta_file_name = "~/<your-workin-directory>/output/AA.fasta")

Interpreting Results from Both Phobius and TMHMM

The usage of TM predictions from both TMHMM and Phobius serves as a powerful tool for evaluating which TM predictions are likely to be accurate biologically. Generally speaking, if TM domains appear in both TMHMM and Phobius predictions, they are likely to be witnessed experimentally as well. If certain domains are only seen on one of the outputs, the likelihood of the domain actually being present is lower. However, as Phobius is able to differentiate between signal peptides and TM domains, any TM regions predicted at the beginning of the sequence by Phobius are likely to be present regardless of TMHMM’s prediction. TMHMM’s TM predictions tend to be more conservative than Phobius’s. Thus, evaluating membrane topology using both methods is an effective way to decide which TM regions are most likely to be present experimentally.

Ranking criteria

To determine which alternatively spliced/mis-spliced isoforms are significantly different from the principal isoforms, surfaltr provides users the option to sort or rank isoform pairs by three different criteria: length, difference in number of TM domains, and a combination metric that takes both length and TM difference into account. The ranking equations are described below:

Choosing and Formatting Isoform Pairs to Display

In addition to ranking the proteins, the plot_isoforms() function can also choose how many proteins to display. Depending on the number supplied to the n_prts input parameter, the user can choose the number of top-ranked pairs to display on the plot. Generally, it is recommended to limit the number of isoform pairs displayed in the plot to under 30 as otherwise the plot can appear cluttered. However, if the user wants to view all isoforms on the same plot, they simply have to enter the total number of genes they entered into the n_prts parameter. As the number of isoform pairs displayed gets larger, the size_txt parameter should be reduced and the space_left parameter should be increased to accommodate more isoform plots. Similarly, if only a few isoforms are plotted, the size_txt parameter should be increased and the space_left parameter decreased. The following table provides helpful guidelines for size parameters that can be used for different values of n_prts:

Plot Interpretation

Although the plots provide an accurate prediction of TM topology for principal and alternative isoforms, there are a few caveats to keep in mind when interpreting the plots. Firstly, on plots created using Phobius’s model, any isoforms that have no TM region and no signal peptide will always be classified as intracellular. Due to the difficulty associated with predicting non-surface protein localization in silico, it is hard to tell if these isoforms are actually intracellular. Thus, these isoforms should not be regarded as intracellular but rather as unknowns. Similarly, on plots created using TMHMM’s model, any isoforms lacking TM regions will be classified as extracellular. Again, it is very difficult to tell computationally if this protein is extracellular and as such, isoforms categorized as fully extracellular should also be regarded as unknowns unless there is experimental information on isoforms from this gene being secreted.

Note: The isoforms pairs seen in the plots produced by surfaltr are not aligned. Therefore, the first amino acid seen in the alternative isoform of one pair may not necessarily correspond to the first amino acid in the principal isoform. In order to facilitate sequence alignment, surfaltr also includes functions that perform multiple sequence alignment of isoform pairs. When used in combination, this can help with the interpretation of isoform pair plots.

From Isoform Pairs from a Single Organism

After identifying genes of interest from an isoform pairs plot, surfaltr offers users an easy way to align amino acid sequences from the same gene through the align_prts() function. Input variables include a vector containing the genes of interest (gene_names), a string containing the path to the file containing the data to be used (data_fie), a boolean stating if the data cotnains amino acids (if_aa), and a string stating which organism the data is from (organism).

The code used to run the align_prts function is shown below:

align_prts(gene_names = c("Crb1", "Adgrl1"), data_file = system.file("extdata", "CRB1.csv", package="surfaltr"), 
           if_aa = TRUE, organism = "mouse")

This function will produce PDF files containing the aligned sequences and FASTA files containing the sequences and transcript IDs for each gene, which will be saved to a folder named ‘output’ in your working directory.

Troubleshooting tips: This function of sufaltr depends on the msa R package which requires LaTeX. Any errors in msa installation pertaining to LaTeX usage and pdf generation can throw errors. If users run into problems please check msa R package isnatllation.

From Isoform Pairs from Multiple Organisms

In addition to single organism alignment, surfaltr also provides users an easy way to align amino acid sequences from the same gene from multiple organisms through the align_org_prts() function. To use align_prts(), the user needs to provide the names of the genes they are interested in examining from both organisms (gene_names), the path to the file cotnaining the human data (hs_data_file), the path to the file containing the mouse data (ms_data_file), and a boolean stating if the data contains amino acid sequences (if_aa).

The code used to run the align_org_prts function is shown below:

align_org_prts(gene_names = c("gene-name"), hs_data_file = "~/path/to/human-dataframe.csv", 
               mm_data_file =  "~/path/to/mouse-dataframe.csv", if_aa = FALSE)

This function will produce PDF files containing the aligned sequences and FASTA files containing the sequences and transcript IDs for each gene, which will be saved to the working directory.

Interpreting Alignment Results and Isoform Pair Plots

library(surfaltr)

sessionInfo()

## R version 4.4.2 (2024-10-31)
## Platform: x86_64-pc-linux-gnu
## Running under: Ubuntu 24.04.1 LTS
## 
## Matrix products: default
## BLAS:   /usr/lib/x86_64-linux-gnu/openblas-pthread/libblas.so.3 
## LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.26.so;  LAPACK version 3.12.0
## 
## locale:
##  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
##  [3] LC_TIME=en_US.UTF-8        LC_COLLATE=C              
##  [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
##  [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
##  [9] LC_ADDRESS=C               LC_TELEPHONE=C            
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       
## 
## time zone: Etc/UTC
## tzcode source: system (glibc)
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] kableExtra_1.4.0 knitr_1.49       surfaltr_1.13.0  testthat_3.2.1.1
## [5] devtools_2.4.5   usethis_3.0.0   
## 
## loaded via a namespace (and not attached):
##   [1] DBI_1.2.3               httr2_1.0.6             remotes_2.5.0          
##   [4] biomaRt_2.63.0          rlang_1.1.4             magrittr_2.0.3         
##   [7] ade4_1.7-22             compiler_4.4.2          RSQLite_2.3.8          
##  [10] systemfonts_1.1.0       png_0.1-8               vctrs_0.6.5            
##  [13] stringr_1.5.1           profvis_0.4.0           pkgconfig_2.0.3        
##  [16] crayon_1.5.3            fastmap_1.2.0           dbplyr_2.5.0           
##  [19] XVector_0.47.0          ellipsis_0.3.2          utf8_1.2.4             
##  [22] promises_1.3.0          rmarkdown_2.29          sessioninfo_1.2.2      
##  [25] tzdb_0.4.0              UCSC.utils_1.3.0        purrr_1.0.2            
##  [28] bit_4.5.0               xfun_0.49               zlibbioc_1.52.0        
##  [31] cachem_1.1.0            seqinr_4.2-36           GenomeInfoDb_1.43.1    
##  [34] jsonlite_1.8.9          progress_1.2.3          blob_1.2.4             
##  [37] later_1.3.2             prettyunits_1.2.0       R6_2.5.1               
##  [40] bslib_0.8.0             stringi_1.8.4           pkgload_1.4.0          
##  [43] brio_1.1.5              jquerylib_0.1.4         Rcpp_1.0.13-1          
##  [46] readr_2.1.5             IRanges_2.41.1          httpuv_1.6.15          
##  [49] tidyselect_1.2.1        rstudioapi_0.17.1       yaml_2.3.10            
##  [52] miniUI_0.1.1.1          curl_6.0.1              pkgbuild_1.4.5         
##  [55] tibble_3.2.1            Biobase_2.67.0          shiny_1.9.1            
##  [58] withr_3.0.2             KEGGREST_1.47.0         evaluate_1.0.1         
##  [61] desc_1.4.3              urlchecker_1.0.1        BiocFileCache_2.15.0   
##  [64] xml2_1.3.6              Biostrings_2.75.1       pillar_1.9.0           
##  [67] filelock_1.0.3          stats4_4.4.2            generics_0.1.3         
##  [70] rprojroot_2.0.4         S4Vectors_0.45.2        hms_1.1.3              
##  [73] ggplot2_3.5.1           munsell_0.5.1           scales_1.3.0           
##  [76] xtable_1.8-4            glue_1.8.0              maketools_1.3.1        
##  [79] tools_4.4.2             sys_3.4.3               buildtools_1.0.0       
##  [82] fs_1.6.5                grid_4.4.2              AnnotationDbi_1.69.0   
##  [85] colorspace_2.1-1        GenomeInfoDbData_1.2.13 msa_1.39.0             
##  [88] cli_3.6.3               rappdirs_0.3.3          protr_1.7-4            
##  [91] fansi_1.0.6             viridisLite_0.4.2       svglite_2.1.3          
##  [94] dplyr_1.1.4             gtable_0.3.6            sass_0.4.9             
##  [97] digest_0.6.37           BiocGenerics_0.53.3     htmlwidgets_1.6.4      
## [100] memoise_2.0.1           htmltools_0.5.8.1       lifecycle_1.0.4        
## [103] httr_1.4.7              mime_0.12               bit64_4.5.2            
## [106] MASS_7.3-61