Dear All,
I am tring to analyse my MeDIP-Seq data in wheat using MEDIPS package. For this purpose, as mentioned in my earlier post, I created .2bit file from .fasta whole genome sequence file and then I created a forge bsgenome package and then build the package using devtools in R. I created medipsets using MEDIPS createSet with extend=300, shift=0, ws=500, paired=T.
The code I used for creating medipset is as follows:
HD0medipnew1=MEDIPS.createSet(file="G:/Lr48/HD0medipnew.bam",BSgenome=BSgenome, paired=TRUE, extend=extend, uniq=uniq, shift=shift, window_size=ws).
Then I tried to generate coupling set using the below code.However, after reaching to chromosome 7, i get the error as follows:
I tried differential coverage with calculating CpG densities and its working by putting MeDIP=F. I would like to know if in case of plants, is it fine to go without calculating CpG densities(coupling set) or not?
I would also like to know if it is required to calculate cpg densities in case of plants or can we go directly for differential coverage without CpG normalisation?
> CS = MEDIPS.couplingVector(pattern = "CG", refObj = HD0medipnew1)
Get genomic sequence pattern positions...
Searching in chr1A ...[ 23590660 ] found.
Searching in chr1B ...[ 27672978 ] found.
Searching in chr1D ...[ 20896261 ] found.
Searching in chr2A ...[ 31176892 ] found.
Searching in chr2B ...[ 32155496 ] found.
Searching in chr2D ...[ 27560821 ] found.
Searching in chr3A ...[ 29825837 ] found.
Searching in chr3B ...[ 33423157 ] found.
Searching in chr3D ...[ 25973362 ] found.
Searching in chr4A ...[ 29541297 ] found.
Searching in chr4B ...[ 27007257 ] found.
Searching in chr4D ...[ 21603933 ] found.
Searching in chr5A ...[ 28266938 ] found.
Searching in chr5B ...[ 101245 ] found.
Searching in chr5D ...[ 175696 ] found.
Searching in chr6A ...[ 300449 ] found.
Searching in chr6B ...[ 436766 ] found.
Searching in chr6D ...[ 528029 ] found.
Searching in chr7A ...Error in .local(con, format, text, ...) : UCSC library operation failed
In addition: Warning message:
In .local(con, format, text, ...) :
twoBitReadSeqFrag in chr7A end (15566174) >= seqSize (15508575)
Could anyone suggest me where I am going wrong?
Your help in this connection will be highly appreciated.
> sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 17134)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.1252 LC_CTYPE=English_United States.1252 LC_MONETARY=English_United States.1252 LC_NUMERIC=C LC_TIME=English_United States.1252
attached base packages:
[1] stats4 parallel stats graphics grDevices utils datasets methods base
other attached packages:
[1] BSgenome.Taestivum.Refseqv1.0_1.0 MEDIPS_1.32.0 Rsamtools_1.32.3 BSgenome_1.48.0 rtracklayer_1.40.6 Biostrings_2.48.0
[7] XVector_0.20.0 GenomicRanges_1.32.7 GenomeInfoDb_1.16.0 IRanges_2.14.12 S4Vectors_0.18.3 BiocGenerics_0.26.0
[13] devtools_1.13.6
loaded via a namespace (and not attached):
[1] Rcpp_0.12.19 BiocInstaller_1.30.0 compiler_3.5.1 BiocManager_1.30.3 git2r_0.23.0 prettyunits_1.0.2 progress_1.2.0 bitops_1.0-6
[9] tools_3.5.1 zlibbioc_1.26.0 biomaRt_2.36.1 digest_0.6.18 bit_1.1-14 preprocessCore_1.42.0 memoise_1.1.0 RSQLite_2.1.1
[17] lattice_0.20-35 pkgconfig_2.0.2 rlang_0.2.2 Matrix_1.2-14 DelayedArray_0.6.6 DBI_1.0.0 GenomeInfoDbData_1.1.0 httr_1.3.1
[25] stringr_1.3.1 withr_2.1.2 hms_0.4.2 gtools_3.8.1 bit64_0.9-7 grid_3.5.1 Biobase_2.40.0 R6_2.3.0
[33] DNAcopy_1.54.0 AnnotationDbi_1.42.1 XML_3.98-1.16 BiocParallel_1.14.2 magrittr_1.5 blob_1.1.1 matrixStats_0.54.0 GenomicAlignments_1.16.0
[41] assertthat_0.2.0 SummarizedExperiment_1.10.1 stringi_1.1.7 RCurl_1.95-4.11 crayon_1.3.4
Regards,
gautam.
Thanks Lukas.
I will try QSEA package as suggested.However, could you point out a possible reason for the error which I encountered during coupling set function?Why the process stopped at reaching to chromosome 7?Do, I need to worry about this or proceed further with my analysis?
Regards,