Graphical options with plot
In this section we show how to plot different graphics with the function plot
. A complete description can be viewed with ?plot.fullsib_scan
.
As an example we will use the data that comes with with the package. First, the steps below are needed; they were described in the fullsibQTL
Tutorial vignette.
library( "fullsibQTL" )
data( "example_QTLfullsib" )
fullsib <- create_fullsib( example_QTLfullsib,
list( LG1_final, LG2_final, LG3_final, LG4_final ),
step = 1,
map.function = "kosambi",
condIndex = 3.5 )
## IM
im <- im_scan( fullsib )
qtls.im <- r2_ls( fullsib,
pheno.col = 1,
lg = c( 1, 1, 2, 2, 2, 3, 3, 4 ),
pos = c( "loc13", "M12", "loc45", "loc77", "loc108", "M33", "loc65", "M52" ) )
## CIM
cofs.fs <- cof_selection( fullsib,
pheno.col = 1,
k=log(300),
n.cofactor = 10)
cim <- cim_scan( cofs.fs )
qtls.cim <- r2_ls( fullsib,
pheno.col = 1,
lg = c( 1, 1, 2, 2, 2, 3, 3, 4 ),
pos = c( "loc13", "M12", "loc45", "loc77", "loc108", "M33", "loc65", "M52" ) )
## Getting im and cim peaks for alpha 0.05.
im.peak1 <- summary( im_perm, alpha = 0.05 )[ 1, 1 ]
im.peak2 <- summary( im_perm, alpha = 0.05 )[ 1, 2 ]
cim.peak1 <- summary( cim_perm, alpha = 0.05 )[ 1, 1 ]
cim.peak2 <- summary( cim_perm, alpha = 0.05 )[ 1, 2 ]
After, to plot a basic graphic:
Customization
The codes below are self-explained. There are a number of modifications in the graphics that can be easily implemented.
## argument lty
## Changes the line type
## There are 6 different line types as follow:
plot( im )
abline( h = c( 1, 2, 3, 4, 5, 6 ), lty = c ( 1, 2, 3, 4, 5, 6 ) )
## argument lwd
### changes the line width. It changes in a relative proportion to
### the default. 2 is twice the default width, 3 is triple of the
### default, and so on.
## You can set any of them with the plot:
plot( im, lty = 3, lwd = 2 )
## argument lg
### selects which LG to plot
## Plotting just LG1
plot( im, lty = 3, lwd = 2, lg = 1 )
## argument label.lg
### labels the LGs
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c( "Chr 1", "Chr 2" ) )
## argument cex.axis
### changes the magnification of axis annotation relative to cex
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c( "Chr 1", "Chr 2" ), cex.axis = 1 )
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c("Chr 1","Chr 2"), cex.axis = 1.5 )
## argument incl.mkr
### includes the marker position there are 3 options:
### `points`, `name`, and, `none`. `none` is default.
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c( "Chr 1", "Chr 2" ), cex.axis = 1.2,
incl.mkr = "points" )
## argument cex.incl
### changes the magnification of incl.mrk annotation
### relative to cex.
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c( "Chr 1", "Chr 2" ), cex.axis = 1.2,
incl.mkr = "name", cex.incl = .75 )
## argument gap
### enlarges the gap between chromosomes (in CM)
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c( "Chr 1", "Chr 2" ), cex.axis = 1.2,
incl.mkr = "name", cex.incl = .75, gap = 50 )
## argument col
### changes the line colors
plot( im, lty = 3, lwd = 2, lg = c( 1, 2 ),
label.lg = c( "Chr 1", "Chr 2" ), cex.axis = 1.2,
incl.mkr = "name", cex.incl = .75, gap = 50, col="red")
You can check more colors options at this cheatsheet or at this color guide
You can also customize the plot with the default R plot.default
options, such as main
, ylab
, xlab
, yaxt
, ylim
, etc.
plot( im, lty=1, lwd=2, lg=c(1,2), cex.axis=1.2,
incl.mkr="name", cex.incl = .9, gap = 50, col="red",
alternate.lgid = TRUE,
## Configuring with R plot options:
main = "IM - Trait - Crop",
xlab = "Linkage Group",
ylab = "LOD score",
yaxt = "n",
ylim = c(0,10))
axis( side = 2, at = seq( 0, 10, 1 ), lwd = 1 )
abline( h = im.peak1, lty = 2 )
abline( h = im.peak2, lty = 3 )
legend( "topleft", lty = c( 1, 2, 3, 0 ), cex = 0.75,
col = c( "red", "black", "black" ),
legend = c( "IM", "1st Peak Threshold", "2nd Peak Threshold" ) )
Graphical comparison between im_scan
and cim_scan
Results from different mapping approaches can be plotted together, with the following codes:
plot( im,cim, lty = c(2,1), lwd = 2, incl.mkr = "name",
cex.incl = 0.7, gap = 50, cex.axis = 0.8,
col = c("blue","red"),
main = "IM - Trait - Crop",
xlab = "Linkage Group",
ylab = "LOD score",
yaxt = "n",
ylim = c( 0, 15 ) )
axis( side = 2, at = seq( 0, 15, 1 ), lwd = 1 )
abline( h = im.peak1, lty = 2, lwd = 0.75 )
abline( h = cim.peak1, lty = 2, lwd = 0.75 )
legend( "topleft", lty = c( 2, 1, 2, 2 ), cex = 0.5,
col = c( "red", "blue", "black", "black" ),
legend = c( "IM", "CIM", "1st Threshold","2nd Threshold") )
Adding QTL locations
The mapped QTLs are the same explained on the fullsibQTL
Tutorial vignette:
print( cim, lg = 1, pos = c( "loc13", "M12" ) )
#> lg pos.cM LOD model
#> loc13 1 13 6.420684 0
#> M12 1 110 5.970667 0
print( cim, lg = 2, pos = c( "loc45", "loc77" ,"loc108" ) )
#> lg pos.cM LOD model
#> loc45 2 45 14.985451 0
#> loc77 2 77 9.115538 0
#> loc108 2 108 10.496445 0
print( cim, lg = 3, pos = c( "M33", "loc65" ) )
#> lg pos.cM LOD model
#> M33 3 20 5.960747 0
#> loc65 3 65 9.367385 0
print( cim, lg = 4, pos = "M52" )
#> lg pos.cM LOD model
#> M52 4 60 11.39792 0
The usage of function plot
to add QTL is not so straightforward as to the one of plot_fullsibQTL
. It is necessary to use function points
as follows. Knowing that the gap
space is 50
, we can add this to the total size of each LG and get the position. The points will be localized at:
points.loc <- c( 13, 110, 45, 77, 108, 20, 65, 60 )
points.loc <- points.loc + c( 0, 0,
140 + 50, 140 + 50, 140 + 50,
140 + 50 + 140 + 50, 140 + 50 + 140 + 50,
140 + 50 + 140 + 50 + 140 + 50 )
And then:
plot( im, cim, lty = c( 2, 1 ), lwd = 2,
incl.mkr = "name", cex.incl = 0.5, gap = 50,
cex.axis = 0.8, col = c( "blue", "red" ),
main = "IM - Trait - Crop",
xlab = "Linkage Group",
ylab = "LOD score",
yaxt = "n",
ylim = c( -2, 15 ) )
axis( side=2, at = seq( 0, 15, 1 ), lwd = 1 )
abline( h = im.peak1, lty = 2, lwd = 0.75 )
abline( h = cim.peak1, lty = 2, lwd = 0.75 )
points( y = rep( 0, 8 ), x = points.loc, pch = 17, col = "dark green" )
legend( "topleft", lty = c( 2, 1, 2, 2, 0 ), cex = 0.75,
pch = c( NA, NA, NA, NA, 17 ),
col = c( "red", "blue", "black", "black", "dark green" ),
legend = c( "IM", "CIM", "1st Threshold", "2nd Threshold", "QTL"),
bty = 'n', y.intersp = .75 )
Session Info
sessionInfo()
#> R version 3.6.3 (2020-02-29)
#> Platform: x86_64-pc-linux-gnu (64-bit)
#> Running under: Ubuntu 19.10
#>
#> Matrix products: default
#> BLAS: /usr/local/lib/R/lib/libRblas.so
#> LAPACK: /usr/local/lib/R/lib/libRlapack.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] fullsibQTL_0.0.9007 onemap_2.1.3 rmdformats_0.3.7
#> [4] knitr_1.28
#>
#> loaded via a namespace (and not attached):
#> [1] colorspace_1.4-1 ellipsis_0.3.1 class_7.3-15
#> [4] rio_0.5.16 htmlTable_1.13.3 base64enc_0.1-3
#> [7] rstudioapi_0.11 mice_3.9.0 codetools_0.2-16
#> [10] splines_3.6.3 doParallel_1.0.15 Formula_1.2-3
#> [13] polynom_1.4-0 jsonlite_1.6.1 MDSMap_1.1
#> [16] broom_0.5.6 cluster_2.1.0 png_0.1-7
#> [19] shiny_1.4.0.2 httr_1.4.1 compiler_3.6.3
#> [22] backports_1.1.7 lazyeval_0.2.2 assertthat_0.2.1
#> [25] Matrix_1.2-18 fastmap_1.0.1 later_1.0.0
#> [28] acepack_1.4.1 htmltools_0.4.0 tools_3.6.3
#> [31] gtable_0.3.0 glue_1.4.1 dplyr_0.8.5
#> [34] Rcpp_1.0.4.6 carData_3.0-3 cellranger_1.1.0
#> [37] vctrs_0.3.0 gdata_2.18.0 nlme_3.1-144
#> [40] iterators_1.0.12 crosstalk_1.1.0.1 xfun_0.13
#> [43] stringr_1.4.0 openxlsx_4.1.5 mime_0.9
#> [46] miniUI_0.1.1.1 lifecycle_0.2.0 weights_1.0.1
#> [49] gtools_3.8.2 princurve_2.1.4 candisc_0.8-3
#> [52] MASS_7.3-51.5 scales_1.1.1 heplots_1.3-5
#> [55] hms_0.5.3 promises_1.1.0 parallel_3.6.3
#> [58] smacof_2.1-0 RColorBrewer_1.1-2 yaml_2.2.1
#> [61] curl_4.3 gridExtra_2.3 ggplot2_3.3.0
#> [64] rpart_4.1-15 reshape_0.8.8 latticeExtra_0.6-29
#> [67] stringi_1.4.6 foreach_1.5.0 plotrix_3.7-8
#> [70] e1071_1.7-3 checkmate_2.0.0 zip_2.0.4
#> [73] manipulateWidget_0.10.1 rlang_0.4.6 pkgconfig_2.0.3
#> [76] rgl_0.100.54 evaluate_0.14 lattice_0.20-38
#> [79] purrr_0.3.4 htmlwidgets_1.5.1 tidyselect_1.1.0
#> [82] plyr_1.8.6 magrittr_1.5 bookdown_0.19
#> [85] R6_2.4.1 generics_0.0.2 nnls_1.4
#> [88] Hmisc_4.4-0 pillar_1.4.4 haven_2.2.0
#> [91] foreign_0.8-75 survival_3.1-8 abind_1.4-5
#> [94] nnet_7.3-12 tibble_3.0.1 crayon_1.3.4
#> [97] car_3.0-7 wordcloud_2.6 plotly_4.9.2.1
#> [100] ellipse_0.4.1 rmarkdown_2.1 jpeg_0.1-8.1
#> [103] grid_3.6.3 readxl_1.3.1 data.table_1.12.8
#> [106] forcats_0.5.0 digest_0.6.25 webshot_0.5.2
#> [109] xtable_1.8-4 tidyr_1.0.3 httpuv_1.5.2
#> [112] munsell_0.5.0 viridisLite_0.3.0