Targeted bisulfite sequencing analysis
Mark Ziemann
2020-12-16
Introduction
In this report, I will take you through an analysis of targeted bisulfite sequencing.
The pipeline used the Biscuit pipeline (https://github.com/zhou-lab/biscuit) run by Larry Croft.
Here I (Mark Ziemann) am reading them into R for analysis with methylKit
Gene regions
GRCh38/hg38 Assembly
chr12:25168373-25168611 POMC1/CASC1
chr12:25161598-25161833 POMC2/CASC1
chr11:78358469-78358619 SOCS/GAB2
chr14:128240139-128241347 LEP window out of range
chr14:24283796-24283936 NPY2/AL096870.8
chr14:24283548-24283843 NPY1/AL096870.8
chr12:31575210-31575467 TNF/DENND5B
CELF2 chr10:10648947-10649946
PPP2R2D chr10:131945290-131946289
LARP4B chr10:864894-865893
GAB2 chr11:78358014-78359013
CASC1 chr12:25161220-25162219
UNKNOWN chr14:128240193-128242193
AL096870 chr14:24283354-24284353
gene_names <- c("CELF2","PPP2R2D","LARP4B","GAB2","CASC1","UNKNOWN","AL096870")
seqnames <- c("10","10","10","11","12","14","14")
start <- c(10648947,131945290,864894,78358014,25161220,128240193,24283354)
end <- c(10649946,131946289,865893,78359013,25162219,128242193,24284353)
genes <- GRanges(seqnames = seqnames, ranges=IRanges(start,end))
genes$gene_names <- gene_namesData preparation
library("GenomicRanges")
library("limma")##
## Attaching package: 'limma'## The following object is masked from 'package:BiocGenerics':
##
## plotMAlibrary("methylKit")
file.list <- list.files(".",pattern="methylkit.tsv")
samplesheet <- read.table("../sample_info/samplesheet_2020-10-15.tsv",sep="\t",header=TRUE)
samplesheet$poolID <- paste("pool",samplesheet$Internal_ID,sep="")
# filter files for those in the sample sheet
file.list <- file.list[which( sapply(strsplit(file.list,"_"),"[[",1) %in% samplesheet$poolID )]
# filter samplesheet for those with files
samplesheet <- samplesheet[which( samplesheet$poolID %in% sapply(strsplit(file.list,"_"),"[[",1) ),]
# reorder the samplesheet
samplenames <- sapply(strsplit(file.list,"_"),"[[",1)
samplesheet <- samplesheet[match(samplenames , samplesheet$poolID),]
samplesheet[,c(25,9)]## poolID GWG
## 20 pool1 20.0
## 11 pool10 16.0
## 2 pool11 12.5
## 3 pool12 13.0
## 15 pool13 17.0
## 12 pool15 16.0
## 9 pool16 15.0
## 24 pool17 23.0
## 4 pool18 13.0
## 6 pool19 13.4
## 21 pool2 20.0
## 5 pool20 13.0
## 16 pool21 17.0
## 17 pool22 17.0
## 19 pool23 18.0
## 22 pool24 20.0
## 18 pool3 17.5
## 14 pool4 17.0
## 1 pool5 12.0
## 13 pool6 16.3
## 10 pool7 15.9
## 23 pool8 22.0
## 7 pool9 13.8# only for debugging
#treatment <- sample(c(0,1),nrow(samplesheet),replace=TRUE)
myobj <- methRead( location=as.list(file.list),
sample.id=as.list(samplesheet$poolID),
pipeline="bismarkCoverage",
header=FALSE,
assembly="hg38",
treatment=samplesheet$GWG,
context="CpG",
mincov = 1 )## Received list of locations.## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.#myobj <- methRead( location=as.list(file.list),
# sample.id=as.list(samplesheet$poolID),
# pipeline="bismarkCoverage",
# header=FALSE,
# assembly="hg38",
# treatment=as.vector(scale(samplesheet$GWG)),
# context="CpG",
# mincov = 1 )
myobjf <- filterByCoverage(myobj,lo.count=5,lo.perc=NULL,
hi.count=NULL,hi.perc=NULL)
meth <- unite(myobj, destrand=FALSE)## uniting...Some basic charts
getCorrelation(meth,plot=TRUE)## pool1 pool10 pool11 pool12 pool13 pool15
## pool1 1.00000000 0.40897691 0.3604163 0.5208813 0.073209045 0.7555899
## pool10 0.40897691 1.00000000 0.4032244 0.3268398 0.333427416 0.2911568
## pool11 0.36041626 0.40322437 1.0000000 0.6198305 0.441029968 0.3959836
## pool12 0.52088130 0.32683980 0.6198305 1.0000000 0.347516924 0.6237518
## pool13 0.07320904 0.33342742 0.4410300 0.3475169 1.000000000 0.1719790
## pool15 0.75558991 0.29115676 0.3959836 0.6237518 0.171978978 1.0000000
## pool16 0.56210831 0.74398122 0.6331427 0.5447539 0.162126625 0.5506619
## pool17 0.24056877 0.08641637 0.2026827 0.3405140 0.013185165 0.2399127
## pool18 0.32994751 0.36595923 0.5699964 0.5163748 0.120825903 0.4150667
## pool19 0.41573786 0.41327550 0.4829392 0.5305814 0.055239235 0.5253873
## pool2 0.57448215 0.45453329 0.5297845 0.6688384 0.387497097 0.7147256
## pool20 0.19772902 0.22135075 0.5461706 0.3934515 0.079351026 0.1993276
## pool21 0.65635905 0.62383074 0.3381257 0.6254302 0.112333272 0.7241945
## pool22 0.73681410 0.51651752 0.3923189 0.5870287 0.012923472 0.8100501
## pool23 0.44187048 0.54924460 0.2591477 0.5465960 -0.005023422 0.3919826
## pool24 0.69156927 0.30469054 0.3531531 0.4984607 0.084739897 0.8193242
## pool3 0.78412954 0.40268211 0.3237503 0.6381204 0.238358652 0.8461336
## pool4 0.48678073 0.64558863 0.5639221 0.4203278 0.296280060 0.4646162
## pool5 0.34081853 0.50920693 0.2776913 0.3833642 -0.045159727 0.2487136
## pool6 0.69136413 0.40606021 0.4728163 0.6844074 0.081626759 0.8300256
## pool7 0.67409015 0.48812677 0.3120739 0.5473428 -0.015485060 0.4882964
## pool8 0.55088164 0.28373719 0.3627691 0.6030270 0.151831711 0.6850121
## pool9 0.61709299 0.45851676 0.6437427 0.5993483 0.157720883 0.6244439
## pool16 pool17 pool18 pool19 pool2 pool20 pool21
## pool1 0.5621083 0.24056877 0.3299475 0.41573786 0.5744822 0.19772902 0.6563591
## pool10 0.7439812 0.08641637 0.3659592 0.41327550 0.4545333 0.22135075 0.6238307
## pool11 0.6331427 0.20268271 0.5699964 0.48293921 0.5297845 0.54617057 0.3381257
## pool12 0.5447539 0.34051397 0.5163748 0.53058143 0.6688384 0.39345148 0.6254302
## pool13 0.1621266 0.01318516 0.1208259 0.05523923 0.3874971 0.07935103 0.1123333
## pool15 0.5506619 0.23991266 0.4150667 0.52538733 0.7147256 0.19932755 0.7241945
## pool16 1.0000000 0.30844212 0.7255991 0.78403559 0.6266400 0.51846370 0.7191023
## pool17 0.3084421 1.00000000 0.2759766 0.38076621 0.4662418 0.15592392 0.1845351
## pool18 0.7255991 0.27597659 1.0000000 0.88023980 0.4408543 0.79937845 0.5208910
## pool19 0.7840356 0.38076621 0.8802398 1.00000000 0.5219483 0.60036477 0.6327935
## pool2 0.6266400 0.46624181 0.4408543 0.52194831 1.0000000 0.15696512 0.7072304
## pool20 0.5184637 0.15592392 0.7993784 0.60036477 0.1569651 1.00000000 0.2158160
## pool21 0.7191023 0.18453508 0.5208910 0.63279354 0.7072304 0.21581599 1.0000000
## pool22 0.7490573 0.39342336 0.5772828 0.71314424 0.7116765 0.33856244 0.8272210
## pool23 0.5132659 0.32166853 0.5688284 0.55404883 0.4108728 0.42269511 0.6716568
## pool24 0.5505269 0.06546836 0.4599826 0.57522855 0.5534148 0.22822519 0.7236396
## pool3 0.5353249 0.40495258 0.3799475 0.51857904 0.7009927 0.17647758 0.7314088
## pool4 0.7060814 0.31860284 0.3765622 0.46144063 0.5938815 0.19949555 0.5572773
## pool5 0.4274494 0.29707411 0.1901085 0.28697061 0.2603578 0.17347030 0.4338619
## pool6 0.7223641 0.42611509 0.6397236 0.72487774 0.7439145 0.36767864 0.7459451
## pool7 0.6539129 0.38792430 0.4745145 0.60053292 0.4774954 0.34963303 0.5871849
## pool8 0.4865401 0.54854075 0.3378517 0.42097719 0.8438003 0.16728084 0.6759903
## pool9 0.8030771 0.32634901 0.8627966 0.82276403 0.5531333 0.68204065 0.6515015
## pool22 pool23 pool24 pool3 pool4 pool5
## pool1 0.73681410 0.441870475 0.69156927 0.7841295 0.4867807 0.34081853
## pool10 0.51651752 0.549244601 0.30469054 0.4026821 0.6455886 0.50920693
## pool11 0.39231889 0.259147745 0.35315306 0.3237503 0.5639221 0.27769127
## pool12 0.58702873 0.546595977 0.49846066 0.6381204 0.4203278 0.38336417
## pool13 0.01292347 -0.005023422 0.08473990 0.2383587 0.2962801 -0.04515973
## pool15 0.81005007 0.391982619 0.81932421 0.8461336 0.4646162 0.24871356
## pool16 0.74905734 0.513265901 0.55052689 0.5353249 0.7060814 0.42744943
## pool17 0.39342336 0.321668533 0.06546836 0.4049526 0.3186028 0.29707411
## pool18 0.57728276 0.568828438 0.45998265 0.3799475 0.3765622 0.19010854
## pool19 0.71314424 0.554048833 0.57522855 0.5185790 0.4614406 0.28697061
## pool2 0.71167652 0.410872839 0.55341481 0.7009927 0.5938815 0.26035776
## pool20 0.33856244 0.422695114 0.22822519 0.1764776 0.1994955 0.17347030
## pool21 0.82722104 0.671656779 0.72363958 0.7314088 0.5572773 0.43386187
## pool22 1.00000000 0.603374640 0.74812847 0.8091796 0.6483278 0.52397384
## pool23 0.60337464 1.000000000 0.40500839 0.5047490 0.3788002 0.62170520
## pool24 0.74812847 0.405008388 1.00000000 0.7191193 0.4305932 0.23364202
## pool3 0.80917958 0.504748984 0.71911928 1.0000000 0.5374803 0.45244656
## pool4 0.64832780 0.378800223 0.43059319 0.5374803 1.0000000 0.67507925
## pool5 0.52397384 0.621705204 0.23364202 0.4524466 0.6750793 1.00000000
## pool6 0.84634063 0.583395029 0.68872435 0.8124908 0.5504218 0.43732364
## pool7 0.64319213 0.603830151 0.44056388 0.6555766 0.5035933 0.52881910
## pool8 0.66934698 0.401527643 0.48387742 0.6817046 0.4844653 0.22174317
## pool9 0.73372716 0.526890915 0.62234147 0.5746223 0.5080477 0.25997571
## pool6 pool7 pool8 pool9
## pool1 0.69136413 0.67409015 0.5508816 0.6170930
## pool10 0.40606021 0.48812677 0.2837372 0.4585168
## pool11 0.47281629 0.31207394 0.3627691 0.6437427
## pool12 0.68440742 0.54734278 0.6030270 0.5993483
## pool13 0.08162676 -0.01548506 0.1518317 0.1577209
## pool15 0.83002561 0.48829638 0.6850121 0.6244439
## pool16 0.72236405 0.65391291 0.4865401 0.8030771
## pool17 0.42611509 0.38792430 0.5485407 0.3263490
## pool18 0.63972357 0.47451446 0.3378517 0.8627966
## pool19 0.72487774 0.60053292 0.4209772 0.8227640
## pool2 0.74391455 0.47749536 0.8438003 0.5531333
## pool20 0.36767864 0.34963303 0.1672808 0.6820406
## pool21 0.74594514 0.58718490 0.6759903 0.6515015
## pool22 0.84634063 0.64319213 0.6693470 0.7337272
## pool23 0.58339503 0.60383015 0.4015276 0.5268909
## pool24 0.68872435 0.44056388 0.4838774 0.6223415
## pool3 0.81249082 0.65557656 0.6817046 0.5746223
## pool4 0.55042184 0.50359334 0.4844653 0.5080477
## pool5 0.43732364 0.52881910 0.2217432 0.2599757
## pool6 1.00000000 0.64225751 0.6565049 0.7299087
## pool7 0.64225751 1.00000000 0.4692330 0.5632958
## pool8 0.65650486 0.46923298 1.0000000 0.4743457
## pool9 0.72990866 0.56329580 0.4743457 1.0000000
clusterSamples(meth, dist="correlation", method="ward", plot=TRUE)## The "ward" method has been renamed to "ward.D"; note new "ward.D2"
##
## Call:
## hclust(d = d, method = HCLUST.METHODS[hclust.method])
##
## Cluster method : ward.D
## Distance : pearson
## Number of objects: 23PCASamples(meth, screeplot=TRUE)
PCASamples(meth)
Statistical analysis
dm <- calculateDiffMeth(meth)## more than two groups detected:
## will calculate methylation difference as the difference of max(x) - min(x),
## where x is vector of mean methylation per group per region,but
## the statistical test will remain the same.## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurreddm## chr start end strand pvalue qvalue meth.diff
## 1 10.10649 10 10649447 * 0.6149363 0.03643923 4.4728021600
## 2 10.10649 10 10649449 * 0.9318240 0.03643923 0.0022105822
## 3 10.10649 10 10649486 * 0.4177847 0.03643923 0.2140592542
## 4 10.10649 10 10649495 * 0.3006589 0.03643923 0.2180915324
## 5 10.10650 10 10649503 * 0.3947909 0.03643923 0.2171792612
## 6 10.10650 10 10649506 * 0.3928871 0.03643923 0.2269705302
## 7 10.10650 10 10649508 * 0.3950802 0.03643923 0.2226717780
## 8 10.13195 10 131945728 * 1.0000000 0.03652597 0.0000000000
## 9 10.13195 10 131945736 * 0.5380315 0.03643923 3.1447636145
## 10 10.13195 10 131945738 * 0.5378048 0.03643923 3.2400279916
## 11 10.13195 10 131945745 * 0.4514108 0.03643923 2.9152317076
## 12 10.13195 10 131945789 * 0.5307665 0.03643923 3.2420638515
## 13 10.86525 10 865252 * 0.2707925 0.03643923 1.8129715655
## 14 10.86525 10 865254 * 0.8710704 0.03643923 0.0005333305
## 15 10.86528 10 865278 * 0.3017122 0.03643923 1.7023072086
## 16 10.86530 10 865296 * 0.2190133 0.03643923 2.0539326800
## 17 10.86535 10 865349 * 0.3645043 0.03643923 1.2858987560
## 18 10.86536 10 865357 * 0.2997304 0.03643923 1.6661993429
## 19 10.86539 10 865391 * 0.3020616 0.03643923 1.6441854414
## 20 10.86539 10 865394 * 0.3023300 0.03643923 1.6755029254
## 21 10.86543 10 865435 * 0.2805800 0.03643923 1.7814073251
## 22 10.86544 10 865441 * 0.2639153 0.03643923 1.7812172558
## 23 10.86547 10 865466 * 0.3014111 0.03643923 1.6672396937
## 24 10.86547 10 865469 * 0.2727660 0.03643923 1.9283532562
## 25 10.86547 10 865473 * 0.2875062 0.03643923 1.7812964803
## 26 10.86548 10 865480 * 0.2830656 0.03643923 1.6553421411
## 27 10.86548 10 865483 * 0.2771621 0.03643923 1.7467021899
## 28 10.86548 10 865485 * 0.2736797 0.03643923 1.7474604603
## 29 11.22056 11 22055661 * 0.8648969 0.03643923 0.1932695968
## 30 11.78359 11 78358500 * 0.8976392 0.03643923 0.1438996924
## 31 11.78359 11 78358506 * 0.8571964 0.03643923 0.1679132618
## 32 11.78359 11 78358514 * 0.8591867 0.03643923 0.1996236017
## 33 11.78359 11 78358539 * 0.9419994 0.03643923 0.0876164617
## 34 11.78359 11 78358556 * 0.9309013 0.03643923 0.1806015528
## 35 11.78359 11 78358581 * 0.9430586 0.03643923 0.1813948024
## 36 11.78359 11 78358587 * 0.8549180 0.03643923 0.1256778565
## 37 12.25162 12 25161616 * 0.7804156 0.03643923 0.8388515567
## 38 12.25162 12 25161619 * 0.6900627 0.03643923 0.6788218592
## 39 12.25162 12 25161623 * 0.7640714 0.03643923 0.6235242476
## 40 12.25162 12 25161628 * 0.7871269 0.03643923 0.5797577259
## 41 12.25162 12 25161637 * 0.8317700 0.03643923 0.6069397799
## 42 12.25162 12 25161653 * 0.8542478 0.03643923 0.5919931727
## 43 12.25162 12 25161661 * 0.9535254 0.03643923 0.5276035376
## 44 12.25162 12 25161671 * 0.8705917 0.03643923 0.5334016644
## 45 12.25162 12 25161686 * 0.8255300 0.03643923 0.5278326073
## 46 12.25162 12 25161689 * 0.7790737 0.03643923 0.4343061451
## 47 12.25162 12 25161699 * 0.7520585 0.03643923 0.3841773201
## 48 12.25162 12 25161716 * 0.7112392 0.03643923 0.6413738534
## 49 12.25162 12 25161720 * 0.7514947 0.03643923 0.6544657794
## 50 12.25162 12 25161728 * 0.7940798 0.03643923 0.5690303930
## 51 12.25162 12 25161742 * 0.7906351 0.03643923 0.6696725697
## 52 12.25162 12 25161764 * 0.7399476 0.03643923 0.6866772419
## 53 12.25162 12 25161766 * 0.7503116 0.03643923 0.6195195181
## 54 12.25162 12 25161785 * 0.7507910 0.03643923 0.6324423789
## 55 12.25162 12 25161792 * 0.7457131 0.03643923 0.6798218760
## 56 12.25162 12 25161799 * 0.7825987 0.03643923 0.6737609680
## 57 12.25162 12 25161804 * 0.7763453 0.03643923 0.7088006106
## 58 12.25162 12 25161811 * 0.7829605 0.03643923 0.6988058086
## 59 12.25162 12 25161813 * 0.8011192 0.03643923 0.6953420601
## 60 12.25162 12 25161817 * 0.8041989 0.03643923 0.6797576615
## 61 12.25162 12 25161825 * 0.8027808 0.03643923 0.7005639629
## 62 14.12824 14 128241074 * 0.7382197 0.03643923 0.4462753175
## 63 14.12824 14 128241076 * 0.6683343 0.03643923 0.4319009388
## 64 14.12824 14 128241078 * 0.8519678 0.03643923 1.2472444324
## 65 14.12824 14 128241090 * 0.7181128 0.03643923 0.4210396392
## 66 14.12824 14 128241092 * 0.6091929 0.03643923 0.2024877061
## 67 14.12824 14 128241095 * 0.6835523 0.03643923 0.2799160252
## 68 14.12824 14 128241108 * 0.8412613 0.03643923 0.0856408792
## 69 14.12824 14 128241116 * 0.7201835 0.03643923 0.4162766558
## 70 14.12824 14 128241118 * 0.9068155 0.03643923 0.1800032728
## 71 14.12824 14 128241132 * 0.8197531 0.03643923 0.9724964751
## 72 14.12824 14 128241151 * 0.8125344 0.03643923 0.4644990046
## 73 14.12824 14 128241155 * 0.8560125 0.03643923 0.1484280124
## 74 14.12824 14 128241160 * 0.9508294 0.03643923 0.4644990046
## 75 14.12824 14 128241163 * 0.6735004 0.03643923 0.2868460593
## 76 14.12824 14 128241178 * 0.7317071 0.03643923 0.2161982372
## 77 14.12824 14 128241183 * 0.8137899 0.03643923 0.8452360344
## 78 14.12824 14 128241193 * 0.6531298 0.03643923 0.1804646966
## 79 14.12824 14 128241204 * 0.8473035 0.03643923 0.4575064866
## 80 14.12824 14 128241207 * 0.7992737 0.03643923 0.1705628574
## 81 14.12824 14 128241216 * 0.8915934 0.03643923 0.4539367401
## 82 14.12824 14 128241227 * 0.6283827 0.03643923 0.7462736009
## 83 14.12824 14 128241240 * 0.7686238 0.03643923 0.3265165836
## 84 14.12824 14 128241245 * 0.6545207 0.03643923 0.3002322551
## 85 14.12824 14 128241259 * 0.6035705 0.03643923 1.0296966193
## 86 14.12824 14 128241277 * 0.8252204 0.03643923 0.4364367063
## 87 14.12824 14 128241286 * 0.5847030 0.03643923 0.9321615823
## 88 14.12824 14 128241291 * 0.8526999 0.03643923 3.8679290554
## 89 14.12824 14 128241297 * 0.7940969 0.03643923 1.2552201218
## 90 14.12824 14 128241316 * 0.8150487 0.03643923 4.2476508947
## 91 14.12824 14 128241319 * 0.5121211 0.03643923 0.9952270061
## 92 14.12824 14 128241321 * 0.6362342 0.03643923 3.2737841200
## 93 14.24284 14 24283574 * 0.6902251 0.03643923 0.2560992047
## 94 14.24284 14 24283578 * 0.8720488 0.03643923 0.0439367311
## 95 14.24284 14 24283596 * 0.9763246 0.03643923 0.3274953723
## 96 14.24284 14 24283618 * 0.9659016 0.03643923 0.3769520732
## 97 14.24284 14 24283628 * 0.6209829 0.03643923 5.2132701422
## 98 14.24284 14 24283641 * 0.8983374 0.03643923 0.9356081453
## 99 14.24284 14 24283672 * 0.9698650 0.03643923 0.9356081453
## 100 14.24284 14 24283696 * 0.9296198 0.03643923 0.3913198150
## 101 14.24284 14 24283698 * 0.6907424 0.03643923 0.9356081453
## 102 14.24284 14 24283710 * 0.5520286 0.03643923 7.6923076923
## 103 14.24284 14 24283713 * 0.8496542 0.03643923 2.7552674230
## 104 14.24284 14 24283722 * 0.8042471 0.03643923 0.4644990046
## 105 14.24284 14 24283724 * 0.9861908 0.03643923 0.9900990099
## 106 14.24284 14 24283731 * 0.9887177 0.03643923 0.3367003367
## 107 14.24284 14 24283776 * 0.7951377 0.03643923 0.9900990099
## 108 14.24284 14 24283820 * 0.9872536 0.03643923 0.0025906065
## 109 14.24284 14 24283854 * 0.9855000 0.03643923 0.0024316417
## 110 14.24284 14 24283859 * 0.9764633 0.03643923 0.0051642033
## 111 14.24284 14 24283894 * 0.9804784 0.03643923 0.0051802537
## 112 14.24284 14 24283905 * 0.9434602 0.03643923 0.0024318899dm$pvalue<0.05## [1] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [13] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [25] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [37] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [49] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [61] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [73] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [85] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [97] FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE FALSE
## [109] FALSE FALSE FALSE FALSEmin(dm$pvalue)## [1] 0.2190133plot(dm$meth.diff,-log10(dm$pvalue) )
hist(dm$meth.diff)
dm.hyper=getMethylDiff(dm,difference=1,qvalue=0.01,type="hyper")## Warning in max(i): no non-missing arguments to max; returning -Infdm.hypo=getMethylDiff(dm,difference=1,qvalue=0.01,type="hypo")## Warning in max(i): no non-missing arguments to max; returning -InfLets run a simple correlation test
z <- samplesheet$GWG
i=0
xx <- apply( percMethylation(meth), 1 , function(x) {
# w <- wilcox.test(z,x)
cor <- cor.test(x,z)
p <- cor$p.value
r <- cor$estimate
i=i+1
HEADER = paste("Site",i," p=",round(p,3)," r=",round(r,3))
plot(z,x,xlab="GWG",ylab="Beta")
mtext(HEADER)
c("p"=p,"r"=r)
})
## Warning in cor(x, y): the standard deviation is zero
xx <- as.data.frame(t(xx))
xx$FDR <- p.adjust(xx[,1],method="fdr")
xx <- data.frame(meth$start,meth$end,xx)
colnames(xx) <- c("seqname","pos","p","r","FDR")
gr <- GRanges(seqnames = xx$seqname, ranges = IRanges(xx$pos) )
gr$r <- xx$r
gr$p <- xx$p
gr$FDR <- xx$FDR
# lets look at each gene
genes[1]$gene_names## [1] "CELF2"ol <- subsetByOverlaps(gr, genes[1])
for ( i in seq_along(genes) ) {
print(genes[i]$gene_names)
ol <- subsetByOverlaps(gr, genes[i])
print( as.data.frame(ol) )
}## [1] "CELF2"
## seqnames start end width strand r p FDR
## 1 10 10649447 10649447 1 * 0.2344185 0.28164639 0.3957310
## 2 10 10649449 10649449 1 * 0.3469867 0.10477161 0.3361232
## 3 10 10649486 10649486 1 * 0.3434430 0.10860695 0.3361232
## 4 10 10649495 10649495 1 * 0.3608112 0.09075757 0.3361232
## 5 10 10649503 10649503 1 * 0.3485078 0.10315628 0.3361232
## 6 10 10649506 10649506 1 * 0.3468705 0.10489573 0.3361232
## 7 10 10649508 10649508 1 * 0.3459315 0.10590307 0.3361232
## [1] "PPP2R2D"
## seqnames start end width strand r p FDR
## 1 10 131945728 131945728 1 * NA NA NA
## 2 10 131945736 131945736 1 * 0.3481118 0.1035750 0.3361232
## 3 10 131945738 131945738 1 * 0.3461035 0.1057181 0.3361232
## 4 10 131945745 131945745 1 * 0.3557843 0.0956816 0.3361232
## 5 10 131945789 131945789 1 * 0.3446872 0.1072488 0.3361232
## [1] "LARP4B"
## seqnames start end width strand r p FDR
## 1 10 865252 865252 1 * -0.3673383 0.08464886 0.3361232
## 2 10 865254 865254 1 * -0.2524397 0.24520076 0.3649784
## 3 10 865278 865278 1 * -0.3626719 0.08898374 0.3361232
## 4 10 865296 865296 1 * -0.4180377 0.04714365 0.3361232
## 5 10 865349 865349 1 * -0.3330116 0.12049309 0.3361232
## 6 10 865357 865357 1 * -0.4064610 0.05427379 0.3361232
## 7 10 865391 865391 1 * -0.3530251 0.09846708 0.3361232
## 8 10 865394 865394 1 * -0.3444896 0.10746372 0.3361232
## 9 10 865435 865435 1 * -0.3651786 0.08663503 0.3361232
## 10 10 865441 865441 1 * -0.3676968 0.08432251 0.3361232
## 11 10 865466 865466 1 * -0.3851019 0.06958379 0.3361232
## 12 10 865469 865469 1 * -0.3926285 0.06385446 0.3361232
## 13 10 865473 865473 1 * -0.3560652 0.09540131 0.3361232
## 14 10 865480 865480 1 * -0.3521908 0.09932107 0.3361232
## 15 10 865483 865483 1 * -0.3666552 0.08527336 0.3361232
## 16 10 865485 865485 1 * -0.3614986 0.09009924 0.3361232
## [1] "GAB2"
## seqnames start end width strand r p FDR
## 1 11 78358500 78358500 1 * 0.2855746 0.1865352 0.3361232
## 2 11 78358506 78358506 1 * 0.2899696 0.1795497 0.3361232
## 3 11 78358514 78358514 1 * 0.2873490 0.1836929 0.3361232
## 4 11 78358539 78358539 1 * 0.2684118 0.2155860 0.3370429
## 5 11 78358556 78358556 1 * 0.2725337 0.2083493 0.3361232
## 6 11 78358581 78358581 1 * 0.2662481 0.2194512 0.3383207
## 7 11 78358587 78358587 1 * 0.2901049 0.1793375 0.3361232
## [1] "CASC1"
## seqnames start end width strand r p FDR
## 1 12 25161616 25161616 1 * 0.2799999 0.1956602 0.3361232
## 2 12 25161619 25161619 1 * 0.2922762 0.1759567 0.3361232
## 3 12 25161623 25161623 1 * 0.2856063 0.1864841 0.3361232
## 4 12 25161628 25161628 1 * 0.2860386 0.1857891 0.3361232
## 5 12 25161637 25161637 1 * 0.2721719 0.2089779 0.3361232
## 6 12 25161653 25161653 1 * 0.2732459 0.2071156 0.3361232
## 7 12 25161661 25161661 1 * 0.2510698 0.2478578 0.3649784
## 8 12 25161671 25161671 1 * 0.2704599 0.2119696 0.3361232
## 9 12 25161686 25161686 1 * 0.2809091 0.1941516 0.3361232
## 10 12 25161689 25161689 1 * 0.2917219 0.1768155 0.3361232
## 11 12 25161699 25161699 1 * 0.3126499 0.1463613 0.3361232
## 12 12 25161716 25161716 1 * 0.2982263 0.1669179 0.3361232
## 13 12 25161720 25161720 1 * 0.2914225 0.1772807 0.3361232
## 14 12 25161728 25161728 1 * 0.2810922 0.1938489 0.3361232
## 15 12 25161742 25161742 1 * 0.2824957 0.1915382 0.3361232
## 16 12 25161764 25161764 1 * 0.2874977 0.1834561 0.3361232
## 17 12 25161766 25161766 1 * 0.2896885 0.1799911 0.3361232
## 18 12 25161785 25161785 1 * 0.2826441 0.1912949 0.3361232
## 19 12 25161792 25161792 1 * 0.2877693 0.1830240 0.3361232
## 20 12 25161799 25161799 1 * 0.2773512 0.2001001 0.3361232
## 21 12 25161804 25161804 1 * 0.2796287 0.1962784 0.3361232
## 22 12 25161811 25161811 1 * 0.2805454 0.1947541 0.3361232
## 23 12 25161813 25161813 1 * 0.2744289 0.2050772 0.3361232
## 24 12 25161817 25161817 1 * 0.2734839 0.2067043 0.3361232
## 25 12 25161825 25161825 1 * 0.2746051 0.2047748 0.3361232
## [1] "UNKNOWN"
## seqnames start end width strand r p FDR
## 1 14 128241074 128241074 1 * -0.29183053 0.17664695 0.3361232
## 2 14 128241076 128241076 1 * -0.22913018 0.29295436 0.4064742
## 3 14 128241078 128241078 1 * 0.12824579 0.55978363 0.6472498
## 4 14 128241090 128241090 1 * -0.15197018 0.48880125 0.5834079
## 5 14 128241092 128241092 1 * -0.50353434 0.01430593 0.3361232
## 6 14 128241095 128241095 1 * -0.36870537 0.08340945 0.3361232
## 7 14 128241108 128241108 1 * -0.25326675 0.24360549 0.3649784
## 8 14 128241116 128241116 1 * -0.10975222 0.61812409 0.7001201
## 9 14 128241118 128241118 1 * -0.28040030 0.19499489 0.3361232
## 10 14 128241132 128241132 1 * 0.33967911 0.11279221 0.3361232
## 11 14 128241151 128241151 1 * 0.03843629 0.86177260 0.8775849
## 12 14 128241155 128241155 1 * -0.31202261 0.14721606 0.3361232
## 13 14 128241160 128241160 1 * 0.24835272 0.25318323 0.3649784
## 14 14 128241163 128241163 1 * -0.48677220 0.01849694 0.3361232
## 15 14 128241178 128241178 1 * -0.31018846 0.14973562 0.3361232
## 16 14 128241183 128241183 1 * 0.08815028 0.68918706 0.7666456
## 17 14 128241193 128241193 1 * -0.41866452 0.04677954 0.3361232
## 18 14 128241204 128241204 1 * 0.04488665 0.83884895 0.8621503
## 19 14 128241207 128241207 1 * -0.21526960 0.32391088 0.4229895
## 20 14 128241216 128241216 1 * 0.06695832 0.76146787 0.8286562
## 21 14 128241227 128241227 1 * -0.19972834 0.36087412 0.4536013
## 22 14 128241240 128241240 1 * -0.29174771 0.17677546 0.3361232
## 23 14 128241245 128241245 1 * -0.43882435 0.03619165 0.3361232
## 24 14 128241259 128241259 1 * -0.19024790 0.38457560 0.4743099
## 25 14 128241277 128241277 1 * 0.07373763 0.73810189 0.8111813
## 26 14 128241286 128241286 1 * -0.17985611 0.41154031 0.5019887
## 27 14 128241291 128241291 1 * 0.21677850 0.32044854 0.4229895
## 28 14 128241297 128241297 1 * 0.05295209 0.81036482 0.8485896
## 29 14 128241316 128241316 1 * 0.20323799 0.35232028 0.4536013
## 30 14 128241319 128241319 1 * -0.19875973 0.36325591 0.4536013
## 31 14 128241321 128241321 1 * 0.08770792 0.69067175 0.7666456
## [1] "AL096870"
## seqnames start end width strand r p FDR
## 1 14 24283574 24283574 1 * -0.171255903 0.43461952 0.5243779
## 2 14 24283578 24283578 1 * -0.049505960 0.82250897 0.8532570
## 3 14 24283596 24283596 1 * 0.021733812 0.92159061 0.9299687
## 4 14 24283618 24283618 1 * -0.003707478 0.98660505 0.9866050
## 5 14 24283628 24283628 1 * 0.327997655 0.12653025 0.3361232
## 6 14 24283641 24283641 1 * 0.237582743 0.27501326 0.3913650
## 7 14 24283672 24283672 1 * 0.133060356 0.54501307 0.6435793
## 8 14 24283696 24283696 1 * -0.055006448 0.80314516 0.8485896
## 9 14 24283698 24283698 1 * 0.289520063 0.18025584 0.3361232
## 10 14 24283710 24283710 1 * 0.220116995 0.31286804 0.4184139
## 11 14 24283713 24283713 1 * 0.249344731 0.25123045 0.3649784
## 12 14 24283722 24283722 1 * 0.128493104 0.55902053 0.6472498
## 13 14 24283724 24283724 1 * -0.119064328 0.58843849 0.6733678
## 14 14 24283731 24283731 1 * -0.058513131 0.79085774 0.8485896
## 15 14 24283776 24283776 1 * -0.428836636 0.04117123 0.3361232
## 16 14 24283820 24283820 1 * 0.198580206 0.36369836 0.4536013
## 17 14 24283854 24283854 1 * 0.282750563 0.19112062 0.3361232
## 18 14 24283859 24283859 1 * 0.224659201 0.30273169 0.4098701
## 19 14 24283894 24283894 1 * 0.224634202 0.30278692 0.4098701
## 20 14 24283905 24283905 1 * 0.335096753 0.11804499 0.3361232Session information
sessionInfo()## R version 4.0.3 (2020-10-10)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.5 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.7.1
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.7.1
##
## 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] parallel stats4 stats graphics grDevices utils datasets
## [8] methods base
##
## other attached packages:
## [1] methylKit_1.16.0 limma_3.46.0 GenomicRanges_1.42.0
## [4] GenomeInfoDb_1.26.0 IRanges_2.24.0 S4Vectors_0.28.0
## [7] BiocGenerics_0.36.0
##
## loaded via a namespace (and not attached):
## [1] mclust_5.4.6 Rcpp_1.0.5
## [3] bdsmatrix_1.3-4 mvtnorm_1.1-1
## [5] lattice_0.20-41 Rsamtools_2.6.0
## [7] Biostrings_2.58.0 gtools_3.8.2
## [9] digest_0.6.27 R6_2.5.0
## [11] plyr_1.8.6 emdbook_1.3.12
## [13] evaluate_0.14 coda_0.19-4
## [15] ggplot2_3.3.2 pillar_1.4.6
## [17] zlibbioc_1.36.0 rlang_0.4.8
## [19] data.table_1.13.2 R.utils_2.10.1
## [21] R.oo_1.24.0 Matrix_1.2-18
## [23] bbmle_1.0.23.1 qvalue_2.22.0
## [25] rmarkdown_2.5 splines_4.0.3
## [27] BiocParallel_1.24.1 stringr_1.4.0
## [29] fastseg_1.36.0 RCurl_1.98-1.2
## [31] munsell_0.5.0 DelayedArray_0.16.0
## [33] rtracklayer_1.50.0 compiler_4.0.3
## [35] numDeriv_2016.8-1.1 xfun_0.19
## [37] pkgconfig_2.0.3 mgcv_1.8-33
## [39] htmltools_0.5.0 SummarizedExperiment_1.20.0
## [41] tidyselect_1.1.0 tibble_3.0.4
## [43] GenomeInfoDbData_1.2.4 matrixStats_0.57.0
## [45] XML_3.99-0.5 crayon_1.3.4
## [47] dplyr_1.0.2 GenomicAlignments_1.26.0
## [49] MASS_7.3-53 bitops_1.0-6
## [51] R.methodsS3_1.8.1 grid_4.0.3
## [53] nlme_3.1-150 gtable_0.3.0
## [55] lifecycle_0.2.0 magrittr_1.5
## [57] scales_1.1.1 KernSmooth_2.23-18
## [59] stringi_1.5.3 XVector_0.30.0
## [61] reshape2_1.4.4 ellipsis_0.3.1
## [63] generics_0.1.0 vctrs_0.3.4
## [65] tools_4.0.3 Biobase_2.50.0
## [67] glue_1.4.2 purrr_0.3.4
## [69] MatrixGenerics_1.2.0 yaml_2.2.1
## [71] colorspace_2.0-0 knitr_1.30