LDHeatMap not displayed

[nickhir]

I am trying to generate a LD heatmap for my plots. When I try out the example, I get the following warning Warning message: Ignoring unknown parameters: max.overlaps, but the heatmap and all other plots get displayed as normal.

eQTpLot(GWAS.df = GWAS.df.example, eQTL.df = eQTL.df.example, gene = "BBS1", 
        gbuild = "hg19", trait = "LDL", tissue =  "Whole_Blood", 
        LD.df =LD.df.example, R2min = 0.25, LDmin = 100)

When I try out my data, I get the same warning message and all plots get displayed normally, the only thing missing is the LD heatmap (the color bar is present). The info says, that it is Generating LDHeatMap with 187 variants

eQTpLot(GWAS.df = gwas_summary, eQTL.df = variant_gene_pairs, Genes.df = gene_df,
        gene = "ENSG00000100350", gbuild = "hg19",  trait = "CLL-PD", 
        tissue =  "all", CollapseMethod = "min",
        sigpvalue_GWAS = 1e-5, range=200, R2min = 0.6,
        LD.df = ld_df, LDmin=7, saveplot = T, getplot = F,
        leadSNP = lead_snp)

Ultimatly the resulting plot looks like this:

Any Idea what might cause this issue?

This is my sessionInfo()

R version 4.0.0 (2020-04-24)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 7 (Core)

Matrix products: default
BLAS/LAPACK: /g/easybuild/x86_64/CentOS/7/haswell/software/OpenBLAS/0.3.9-GCC-9.3.0/lib/libopenblas_haswellp-r0.3.9.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [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       

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] eQTpLot_0.0.0.9000 forcats_0.5.0      stringr_1.4.0      dplyr_1.0.7       
 [5] purrr_0.3.4        readr_1.3.1        tidyr_1.0.2        tibble_3.0.1      
 [9] ggplot2_3.3.0      tidyverse_1.3.0    data.table_1.12.8 

loaded via a namespace (and not attached):
  [1] readxl_1.3.1                backports_1.1.6            
  [3] Hmisc_4.4-0                 BiocFileCache_1.14.0       
  [5] lazyeval_0.2.2              splines_4.0.0              
  [7] BiocParallel_1.24.1         GenomeInfoDb_1.26.7        
  [9] digest_0.6.25               ensembldb_2.14.1           
 [11] htmltools_0.4.0             fansi_0.4.1                
 [13] magrittr_1.5                checkmate_2.0.0            
 [15] memoise_1.1.0               BSgenome_1.58.0            
 [17] cluster_2.1.0               openxlsx_4.1.4             
 [19] Biostrings_2.58.0           modelr_0.1.6               
 [21] matrixStats_0.59.0          R.utils_2.9.2              
 [23] askpass_1.1                 prettyunits_1.1.1          
 [25] jpeg_0.1-8.1                colorspace_1.4-1           
 [27] blob_1.2.1                  rvest_0.3.5                
 [29] rappdirs_0.3.1              ggrepel_0.8.2              
 [31] haven_2.2.0                 xfun_0.24                  
 [33] crayon_1.3.4                RCurl_1.98-1.2             
 [35] jsonlite_1.6.1              survival_3.1-12            
 [37] VariantAnnotation_1.36.0    glue_1.4.1.9000            
 [39] gtable_0.3.0                zlibbioc_1.36.0            
 [41] XVector_0.30.0              DelayedArray_0.16.3        
 [43] car_3.0-7                   BiocGenerics_0.36.1        
 [45] abind_1.4-5                 scales_1.1.0               
 [47] DBI_1.1.0                   rstatix_0.5.0              
 [49] Rcpp_1.0.4.6                viridisLite_0.3.0          
 [51] progress_1.2.2              htmlTable_1.13.3           
 [53] gridGraphics_0.5-0          foreign_0.8-79             
 [55] bit_1.1-15.2                Formula_1.2-3              
 [57] stats4_4.0.0                htmlwidgets_1.5.1          
 [59] httr_1.4.1                  RColorBrewer_1.1-2         
 [61] LDheatmap_1.0-4             acepack_1.4.1              
 [63] ellipsis_0.3.2              pkgconfig_2.0.3            
 [65] XML_3.99-0.3                R.methodsS3_1.8.0          
 [67] farver_2.0.3                Gviz_1.32.0                
 [69] nnet_7.3-14                 dbplyr_1.4.3               
 [71] utf8_1.1.4                  labeling_0.3               
 [73] ggplotify_0.0.5             tidyselect_1.1.1           
 [75] rlang_0.4.11                AnnotationDbi_1.52.0       
 [77] munsell_0.5.0               cellranger_1.1.0           
 [79] tools_4.0.0                 cli_3.0.0                  
 [81] generics_0.0.2              RSQLite_2.2.0              
 [83] broom_0.5.6                 knitr_1.28                 
 [85] bit64_0.9-7                 fs_1.4.1                   
 [87] zip_2.0.4                   AnnotationFilter_1.14.0    
 [89] nlme_3.1-147                R.oo_1.23.0                
 [91] xml2_1.3.2                  biomaRt_2.46.3             
 [93] compiler_4.0.0              rstudioapi_0.11            
 [95] curl_4.3                    png_0.1-7                  
 [97] ggsignif_0.6.0              reprex_0.3.0               
 [99] stringi_1.4.6               GenomicFeatures_1.42.3     
[101] lattice_0.20-41             ProtGenerics_1.22.0        
[103] Matrix_1.2-18               vctrs_0.3.8                
[105] pillar_1.6.1                lifecycle_1.0.0            
[107] BiocManager_1.30.16         bitops_1.0-6               
[109] patchwork_1.0.0             rtracklayer_1.50.0         
[111] GenomicRanges_1.42.0        R6_2.4.1                   
[113] latticeExtra_0.6-29         rio_0.5.16                 
[115] gridExtra_2.3               IRanges_2.24.1             
[117] dichromat_2.0-0             assertthat_0.2.1           
[119] SummarizedExperiment_1.20.0 openssl_1.4.1              
[121] withr_2.2.0                 GenomicAlignments_1.26.0   
[123] Rsamtools_2.6.0             S4Vectors_0.28.1           
[125] GenomeInfoDbData_1.2.4      mgcv_1.8-31                
[127] parallel_4.0.0              hms_0.5.3                  
[129] grid_4.0.0                  rpart_4.1-15               
[131] rvcheck_0.1.8               carData_3.0-3              
[133] MatrixGenerics_1.2.1        ggpubr_0.3.0               
[135] biovizBase_1.38.0           Biobase_2.50.0             
[137] lubridate_1.7.8             base64enc_0.1-3 

[anastasia-lucas]

Hi, I have not come across this issue before. Would you be able to provide a reproducible example, i.e. either a small subset of your data or some dummy data for GWAS.df, eQTL.df, Genes.df, LD.df, and lead_snp that result in the same issue?

[nickhir]

Alright, I have uploaded the data frames that you requested. You can find them here: https://we.tl/t-qaj0qZZxGm

This is the code I use to generate the plot afterwards:

library(eQTpLot)

load("Genes.df.rds")
load("LD.df.rds")
load("eQTL.df.rds")
load("GWAS.df.rds")
lead_snp <- "rs112378630"
options(bitmapType="cairo-png") # necessary, otherwise error unable to open connection to X11 display
eQTpLot(GWAS.df = subset_gwas_summary, eQTL.df = eQTL.df, Genes.df = gene_df,
        gene = "ENSG00000100350", gbuild = "hg19",  trait = "X", 
        tissue =  "all", CollapseMethod = "min",
        sigpvalue_GWAS = 1e-5, range=200, R2min = 0.6,
        LD.df = LD.df, LDmin=7, saveplot =T, getplot = F,
        leadSNP = lead_snp)

[anastasia-lucas]

Perfect, thank you @nickhir I was able to reproduce the error and am working on finding the issue.

[anastasia-lucas]

Hi @nickhir

There was a bug in how the script was performing the LD filtering. I have tried your example and it was able to produce the correct output.

The repo should be now updated with the fix and you can get the latest version by reinstalling the package by using devtools::install_github('RitchieLab/eQTpLot')

I will close the issue, but please feel free to reopen if you still have another problem.

[nickhir]

Thank you very much for taking care of this issue! The plot now gets generated. However, I am still a little confused by the output and not quite sure if it gets displayed as intended. If I modify the LD.df artificially to only plot a handful of SNPs the LD heatmap does not look as expected. This is the code I used:

library(eQTpLot)

load("Genes.df.rds")
load("LD.df.rds")
load("eQTL.df.rds")
load("GWAS.df.rds")
lead_snp <- "rs112378630"

# select a subset for demonstration purpose 
LD_subset <- LD.df %>%
  filter(R2 > 0.9) %>%
  filter(BP_A > 36783555) %>%
  filter(BP_A < 36800000)

options(bitmapType="cairo-png") # necessary, otherwise error unable to open connection to X11 display
eQTpLot(GWAS.df = subset_gwas_summary, eQTL.df = eQTL.df, Genes.df = gene_df,
        gene = "ENSG00000100350", gbuild = "hg19",  trait = "X", 
        tissue =  "all", CollapseMethod = "min",
        sigpvalue_GWAS = 1e-5, range=200, R2min = 0.0,
        LD.df = LD_subset, LDmin=0, saveplot =F, getplot = T,
        leadSNP = lead_snp)

The resulting imagae looks something like this

If I understand the display correctly, the SNPs do not align with their actual positions.

[anastasia-lucas]

@nickhir I think the latest update https://github.com/RitchieLab/eQTpLot/commit/b46ff6772d2ecd6258ef98b6b81fab019c536049 should resolve this issue. I will wait to close until you can confirm it's working now.

[nickhir]

Unfortunately it still seems like the alignment is not right. It is closer to the real position but still off by a bit. Also, the heatmap doesn't get displayed anymore, similar to the initial problem:

As you can see, the lines do not line up with the lead SNP. The output of LD_subset is:

  CHR_A     BP_A       SNP_A CHR_B     BP_B      SNP_B       R2
1    22 36783638 rs112378630    22 36786966  rs6000269 1.000000
2    22 36788357   rs7290251    22 36790715 rs62228884 0.975001
3    22 36788357   rs7290251    22 36791833  rs7364216 0.975001
4    22 36789003  rs62228881    22 36790234 rs59177520 1.000000

This is the code I used to generate the plot:

devtools::install_github('RitchieLab/eQTpLot')
library(eQTpLot)

load("Genes.df.rds")
load("LD.df.rds")
load("eQTL.df.rds")
load("GWAS.df.rds")
lead_snp <- "rs112378630"

# select a subset for demonstration purpose 
LD_subset <- LD.df %>%
  filter(R2 > 0.9) %>%
  filter(BP_A > 36783555) %>%
  filter(BP_A < 36790000)

options(bitmapType="cairo-png") # necessary, otherwise error unable to open connection to X11 display
eQTpLot(GWAS.df = subset_gwas_summary, eQTL.df = eQTL.df, Genes.df = gene_df,
        gene = "ENSG00000100350", gbuild = "hg19",  trait = "X", 
        tissue =  "all", CollapseMethod = "min",
        sigpvalue_GWAS = 1e-5, range=200, R2min = 0.0,
        LD.df = LD_subset, LDmin=0, saveplot =F, getplot = T,
        leadSNP = lead_snp)

This is the package description.

packageDescription('eQTpLot')
License: GPL-3
Encoding: UTF-8
LazyData: true
Roxygen: list(markdown = TRUE)
RoxygenNote: 7.1.1
Imports: dplyr, patchwork, gridExtra, Gviz, ggnewscale, GenomicRanges, biomaRt, ggpubr, ggplotify, LDheatmap, viridisLite
Depends: ggplot2 (>= 3.3.0), R (>= 3.5.0)
RemoteType: github
RemoteHost: api.github.com
RemoteRepo: eQTpLot
RemoteUsername: RitchieLab
RemoteRef: master
RemoteSha: b46ff6772d2ecd6258ef98b6b81fab019c536049
GithubRepo: eQTpLot
GithubUsername: RitchieLab
GithubRef: master
GithubSHA1: b46ff6772d2ecd6258ef98b6b81fab019c536049
NeedsCompilation: no
Packaged: 2021-08-18 07:31:38 UTC; 
Author: Theodore Drivas [aut], Anastasia Lucas [cre]
Maintainer: Anastasia Lucas <anastasia.lucas.bioinfo@gmail.com>
Built: R 4.0.0; ; 2021-08-18 07:31:40 UTC; unix

[anastasia-lucas]

@nickhir I have pushed an update https://github.com/RitchieLab/eQTpLot/commit/dd8bc290511602908b4b75f923d802fe6ddc1517 that should hopefully improve the alignment.

Unfortunately I was not able to replicate your issue of the heatmap not appearing. This is the resulting figure I get from your code and data (note I did not include the bitmap options on my system), across a few versions of R.

[nickhir]

Thank you for your continuous support. It is very appreciated!