Contents

1 Introduction

This R package provides methods for genetic finemapping in inbred mice by taking advantage of their very high homozygosity rate (>95%).

For one ore more chromosomal regions (GRCm38), method finemap extracts homozygous SNVs for which the allele differs between two sets of strains (e.g. case vs controls) and outputs respective causal SNV/gene candidates.

2 Installation

if(!requireNamespace("BiocManager", quietly = TRUE))
    install.packages("BiocManager")
BiocManager::install("MouseFM")

3 Loading package

library(MouseFM)

4 Example function calls

Available mouse strains

avail_strains()
#>              id        strain
#> 1  129P2_OlaHsd  129P2/OlaHsd
#> 2   129S1_SvImJ   129S1/SvImJ
#> 3  129S5SvEvBrd 129S5/SvEvBrd
#> 4           A_J           A/J
#> 5         AKR_J         AKR/J
#> 6       BALB_cJ       BALB/cJ
#> 7          BTBR          BTBR
#> 8       BUB_BnJ       BUB/BnJ
#> 9       C3H_HeH       C3H/HeH
#> 10      C3H_HeJ       C3H/HeJ
#> 11    C57BL_10J     C57BL/10J
#> 12     C57BL_6J      C57BL/6J
#> 13    C57BL_6NJ     C57BL/6NJ
#> 14    C57BR_cdJ     C57BR/cdJ
#> 15       C57L_J        C57L/J
#> 16        C58_J         C58/J
#> 17     CAST_EiJ      CAST/EiJ
#> 18        CBA_J         CBA/J
#> 19       DBA_1J        DBA/1J
#> 20       DBA_2J        DBA/2J
#> 21       FVB_NJ        FVB/NJ
#> 22        I_LnJ         I/LnJ
#> 23       KK_HiJ        KK/HiJ
#> 24    LEWES_EiJ     LEWES/EiJ
#> 25         LP_J          LP/J
#> 26     MOLF_EiJ      MOLF/EiJ
#> 27   NOD_ShiLtJ    NOD/ShiLtJ
#> 28     NZB_B1NJ      NZB/B1NJ
#> 29    NZO_HlLtJ     NZO/HlLtJ
#> 30     NZW_LacJ      NZW/LacJ
#> 31      PWK_PhJ       PWK/PhJ
#> 32         RF_J          RF/J
#> 33      SEA_GnJ       SEA/GnJ
#> 34    SPRET_EiJ     SPRET/EiJ
#> 35        ST_bJ         ST/bJ
#> 36      WSB_EiJ       WSB/EiJ
#> 37  ZALENDE_EiJ   ZALENDE/EiJ

Call finemap to finemap a specific region

res = finemap(chr="chr7",
              strain1=c("C57BL_6J","C57L_J","CBA_J","NZB_B1NJ"),
              strain2=c("C3H_HEJ","MOLF_EiJ","NZW_LacJ","WSB_EiJ","SPRET_EiJ"),
              impact=c("HIGH", "MODERATE", "LOW"))
#> Query chr7

res[1:10,]
#>    chr      pos        rsid ref alt
#> 1    7 45666192  rs51324364   A   G
#> 2    7 45853238  rs47469186   T   C
#> 3    7 45858570  rs47348864   A   C
#> 4    7 45977282  rs32753716   A   G
#> 5    7 45996764  rs32757904   T   C
#> 6    7 45996772  rs51886013   A   G
#> 7    7 45998716  rs32753986   A   G
#> 8    7 46029114  rs46389823   A   G
#> 9    7 46068710  rs32761224   A   C
#> 10   7 46081400 rs108318096   T   C
#>                                                                                             consequences
#> 1  non_coding_transcript_exon_variant,non_coding_transcript_variant,intron_variant,splice_region_variant
#> 2                                                                                     synonymous_variant
#> 3                                                                   intron_variant,splice_region_variant
#> 4    upstream_gene_variant,non_coding_transcript_exon_variant,downstream_gene_variant,synonymous_variant
#> 5                                                    non_coding_transcript_exon_variant,missense_variant
#> 6                                                  non_coding_transcript_exon_variant,synonymous_variant
#> 7                                                  non_coding_transcript_exon_variant,synonymous_variant
#> 8                                                               upstream_gene_variant,synonymous_variant
#> 9                                                                                       missense_variant
#> 10                                                              upstream_gene_variant,synonymous_variant
#>    C57BL_6J C3H_HeJ C57L_J CBA_J MOLF_EiJ NZB_B1NJ NZW_LacJ SPRET_EiJ WSB_EiJ
#> 1         0       1      0     0        1        0        1         1       1
#> 2         0       1      0     0        1        0        1         1       1
#> 3         0       1      0     0        1        0        1         1       1
#> 4         0       1      0     0        1        0        1         1       1
#> 5         0       1      0     0        1        0        1         1       1
#> 6         0       1      0     0        1        0        1         1       1
#> 7         0       1      0     0        1        0        1         1       1
#> 8         0       1      0     0        1        0        1         1       1
#> 9         0       1      0     0        1        0        1         1       1
#> 10        0       1      0     0        1        0        1         1       1

View meta information

comment(res)
#> [1] "#Alleles of strain C57BL_6J represent the reference (ref) alleles"
#> [2] "#reference_version=GRCm38"

Annotate with consequences (Not recommended for large queries!)

cons = annotate_consequences(res, "mouse")

Annotate with genes

genes = annotate_mouse_genes(res, flanking = 50000)

5 Output of Session Info

The output of sessionInfo() on the system on which this document was compiled:

sessionInfo()
#> R version 4.2.1 (2022-06-23)
#> Platform: x86_64-pc-linux-gnu (64-bit)
#> Running under: Ubuntu 20.04.5 LTS
#> 
#> Matrix products: default
#> BLAS:   /home/biocbuild/bbs-3.16-bioc/R/lib/libRblas.so
#> LAPACK: /home/biocbuild/bbs-3.16-bioc/R/lib/libRlapack.so
#> 
#> locale:
#>  [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
#>  [3] LC_TIME=en_GB              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       
#> 
#> attached base packages:
#> [1] stats     graphics  grDevices utils     datasets  methods   base     
#> 
#> other attached packages:
#> [1] MouseFM_1.8.0    BiocStyle_2.26.0
#> 
#> loaded via a namespace (and not attached):
#>  [1] Biobase_2.58.0         httr_1.4.4             tidyr_1.2.1           
#>  [4] sass_0.4.2             bit64_4.0.5            jsonlite_1.8.3        
#>  [7] gtools_3.9.3           bslib_0.4.0            assertthat_0.2.1      
#> [10] BiocManager_1.30.19    stats4_4.2.1           BiocFileCache_2.6.0   
#> [13] blob_1.2.3             GenomeInfoDbData_1.2.9 yaml_2.3.6            
#> [16] progress_1.2.2         pillar_1.8.1           RSQLite_2.2.18        
#> [19] rlist_0.4.6.2          glue_1.6.2             digest_0.6.30         
#> [22] GenomicRanges_1.50.0   XVector_0.38.0         colorspace_2.0-3      
#> [25] plyr_1.8.7             htmltools_0.5.3        XML_3.99-0.12         
#> [28] pkgconfig_2.0.3        biomaRt_2.54.0         bookdown_0.29         
#> [31] zlibbioc_1.44.0        purrr_0.3.5            scales_1.2.1          
#> [34] tibble_3.1.8           KEGGREST_1.38.0        generics_0.1.3        
#> [37] IRanges_2.32.0         ggplot2_3.3.6          ellipsis_0.3.2        
#> [40] cachem_1.0.6           BiocGenerics_0.44.0    cli_3.4.1             
#> [43] magrittr_2.0.3         crayon_1.5.2           memoise_2.0.1         
#> [46] evaluate_0.17          fansi_1.0.3            xml2_1.3.3            
#> [49] tools_4.2.1            data.table_1.14.4      prettyunits_1.1.1     
#> [52] hms_1.1.2              lifecycle_1.0.3        stringr_1.4.1         
#> [55] S4Vectors_0.36.0       munsell_0.5.0          AnnotationDbi_1.60.0  
#> [58] Biostrings_2.66.0      compiler_4.2.1         jquerylib_0.1.4       
#> [61] GenomeInfoDb_1.34.0    rlang_1.0.6            grid_4.2.1            
#> [64] RCurl_1.98-1.9         rappdirs_0.3.3         bitops_1.0-7          
#> [67] rmarkdown_2.17         gtable_0.3.1           DBI_1.1.3             
#> [70] curl_4.3.3             reshape2_1.4.4         R6_2.5.1              
#> [73] knitr_1.40             dplyr_1.0.10           fastmap_1.1.0         
#> [76] bit_4.0.4              utf8_1.2.2             filelock_1.0.2        
#> [79] stringi_1.7.8          Rcpp_1.0.9             vctrs_0.5.0           
#> [82] png_0.1-7              dbplyr_2.2.1           tidyselect_1.2.0      
#> [85] xfun_0.34