This code is available at https://github.com/aaronsk7/guppy-methylation

This script performs differential methylation analysis using MethylKit.

library("R.utils")
## Loading required package: R.oo
## Loading required package: R.methodsS3
## R.methodsS3 v1.8.1 (2020-08-26 16:20:06 UTC) successfully loaded. See ?R.methodsS3 for help.
## R.oo v1.24.0 (2020-08-26 16:11:58 UTC) successfully loaded. See ?R.oo for help.
## 
## Attaching package: 'R.oo'
## The following object is masked from 'package:R.methodsS3':
## 
##     throw
## The following objects are masked from 'package:methods':
## 
##     getClasses, getMethods
## The following objects are masked from 'package:base':
## 
##     attach, detach, load, save
## R.utils v2.11.0 (2021-09-26 08:30:02 UTC) successfully loaded. See ?R.utils for help.
## 
## Attaching package: 'R.utils'
## The following object is masked from 'package:utils':
## 
##     timestamp
## The following objects are masked from 'package:base':
## 
##     cat, commandArgs, getOption, inherits, isOpen, nullfile, parse,
##     warnings
library("parallel")

library("reshape2")
# library("kableExtra") having some trouble installing dependancies on bio2
library("dplyr")
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library("RColorBrewer")
library("GenomicRanges")
## Loading required package: stats4
## Loading required package: BiocGenerics
## 
## Attaching package: 'BiocGenerics'
## The following objects are masked from 'package:dplyr':
## 
##     combine, intersect, setdiff, union
## The following objects are masked from 'package:parallel':
## 
##     clusterApply, clusterApplyLB, clusterCall, clusterEvalQ,
##     clusterExport, clusterMap, parApply, parCapply, parLapply,
##     parLapplyLB, parRapply, parSapply, parSapplyLB
## The following objects are masked from 'package:stats':
## 
##     IQR, mad, sd, var, xtabs
## The following objects are masked from 'package:base':
## 
##     anyDuplicated, append, as.data.frame, basename, cbind, colnames,
##     dirname, do.call, duplicated, eval, evalq, Filter, Find, get, grep,
##     grepl, intersect, is.unsorted, lapply, Map, mapply, match, mget,
##     order, paste, pmax, pmax.int, pmin, pmin.int, Position, rank,
##     rbind, Reduce, rownames, sapply, setdiff, sort, table, tapply,
##     union, unique, unsplit, which, which.max, which.min
## Loading required package: S4Vectors
## 
## Attaching package: 'S4Vectors'
## The following objects are masked from 'package:dplyr':
## 
##     first, rename
## The following object is masked from 'package:base':
## 
##     expand.grid
## Loading required package: IRanges
## 
## Attaching package: 'IRanges'
## The following objects are masked from 'package:dplyr':
## 
##     collapse, desc, slice
## The following object is masked from 'package:R.oo':
## 
##     trim
## Loading required package: GenomeInfoDb
library("limma")
## 
## Attaching package: 'limma'
## The following object is masked from 'package:BiocGenerics':
## 
##     plotMA
library("methylKit")
## 
## Attaching package: 'methylKit'
## The following object is masked from 'package:dplyr':
## 
##     select
library("gplots")
## 
## Attaching package: 'gplots'
## The following object is masked from 'package:IRanges':
## 
##     space
## The following object is masked from 'package:S4Vectors':
## 
##     space
## The following object is masked from 'package:stats':
## 
##     lowess
library("seqinr")
## 
## Attaching package: 'seqinr'
## The following object is masked from 'package:limma':
## 
##     zscore
## The following object is masked from 'package:dplyr':
## 
##     count
## The following object is masked from 'package:R.oo':
## 
##     getName

Import biscuit data with a function

Define and test the function

vcfz <- "meth_data/test-R1.bam.vcf.gz.cg.bed.gz"

biscuit2methylkit <- function(vcfz) {
vcf <- gsub(".gz$","",vcfz)
gunzip(vcfz)
tdata <- read.table(vcf)
gzip(vcf)
chrBase <- paste(tdata$V1,tdata$V2,sep=".")
chr <- tdata$V1
base <- tdata$V2
strand="F"
coverage <- tdata$V5
freqC <- tdata$V4 *100 
freqT <- 100-freqC 
mk <- data.frame(  chrBase , chr , base, strand,  coverage, freqC, freqT )
mk
mkname <- paste(vcf,".mk",sep="")
mkname <- gsub(".bam.vcf.gz.cg.bed","",mkname)
write.table(mk,file=mkname,quote=FALSE,sep="\t",row.names = FALSE)
}


biscuit2methylkit(vcfz)
getwd()
## [1] "/home/ziemannm/guppy_methylseq/main_data1"
myobj=methRead("meth_data/test-R1.mk",
           sample.id="test1",
           assembly="hg18",
           treatment=1,
           context="CpG",
           mincov = 10 )
## Received single location.
## Reading file.
myobj
##            chr  start    end strand coverage numCs numTs
## 1   KK215283.1    113    113      +      126     3   123
## 2   KK215283.1   2071   2071      +       10     8     2
## 3   KK215283.1   2111   2111      +       13    11     2
## 4   KK215283.1   2139   2139      +       13    13     0
## 5   KK215283.1   4814   4814      +       10     8     2
## 6   KK215283.1   5196   5196      +       10    10     0
## 7   KK215283.1   5241   5241      +       10     7     3
## 8   KK215283.1   5882   5882      +       10    10     0
## 9   KK215283.1   5885   5885      +       10    10     0
## 10  KK215283.1   7324   7324      +       10     9     1
## 11  KK215283.1   8493   8493      +       10     0    10
## 12  KK215283.1   8605   8605      +       10     9     1
## 13  KK215283.1   8617   8617      +       11    10     1
## 14  KK215283.1   8814   8814      +       10     9     1
## 15  KK215283.1  11557  11557      +       15    13     2
## 16  KK215283.1  11638  11638      +       12    12     0
## 17  KK215283.1  12468  12468      +       12    12     0
## 18  KK215283.1  12497  12497      +       10    10     0
## 19  KK215283.1  12506  12506      +       10     9     1
## 20  KK215283.1  12531  12531      +       11    10     1
## 21  KK215283.1  12560  12560      +       11     9     2
## 22  KK215283.1  12575  12575      +       12    10     2
## 23  KK215283.1  12580  12580      +       12     7     5
## 24  KK215283.1  12586  12586      +       12     7     5
## 25  KK215283.1  15071  15071      +       12     8     4
## 26  KK215283.1  16099  16099      +       10    10     0
## 27  KK215283.1  16120  16120      +       11    10     1
## 28  KK215283.1  16149  16149      +       11    11     0
## 29  KK215283.1  16175  16175      +       13    10     3
## 30  KK215283.1  16188  16188      +       10     7     3
## 31  KK215283.1  16647  16647      +       10     8     2
## 32  KK215283.1  16675  16675      +       10     7     3
## 33  KK215283.1  16980  16980      +       10    10     0
## 34  KK215283.1  17746  17746      +       11     9     2
## 35  KK215283.1  17771  17771      +       10     8     2
## 36  KK215283.1  21884  21884      +       10     6     4
## 37  KK215283.1  22800  22800      +       12    12     0
## 38  KK215283.1  26664  26664      +       10     7     3
## 39  KK215283.1  27193  27193      +       11    10     1
## 40  KK215283.1  27205  27205      +       10     9     1
## 41  KK215283.1  27224  27224      +       10     9     1
## 42  KK215283.1  27234  27234      +       10     8     2
## 43  KK215283.1  27978  27978      +       10     9     1
## 44  KK215283.1  29546  29546      +       13    13     0
## 45  KK215283.1  29628  29628      +       11     6     5
## 46  KK215283.1  30232  30232      +       13    12     1
## 47  KK215283.1  30277  30277      +       10     9     1
## 48  KK215283.1  30776  30776      +       10     4     6
## 49  KK215283.1  30909  30909      +       11     9     2
## 50  KK215283.1  31777  31777      +       10    10     0
## 51  KK215283.1  31779  31779      +       10    10     0
## 52  KK215283.1  31950  31950      +       11    11     0
## 53  KK215283.1  31961  31961      +       11    11     0
## 54  KK215283.1  31970  31970      +       10     7     3
## 55  KK215283.1  32001  32001      +       10    10     0
## 56  KK215283.1  32004  32004      +       10    10     0
## 57  KK215283.1  32030  32030      +       10     9     1
## 58  KK215283.1  32038  32038      +       11    11     0
## 59  KK215283.1  32042  32042      +       11    11     0
## 60  KK215283.1  32092  32092      +       10     8     2
## 61  KK215283.1  32101  32101      +       11    10     1
## 62  KK215283.1  32703  32703      +       10     7     3
## 63  KK215283.1  32958  32958      +       10     1     9
## 64  KK215283.1  32963  32963      +       10     3     7
## 65  KK215283.1  32970  32970      +       11     2     9
## 66  KK215283.1  32976  32976      +       10     1     9
## 67  KK215283.1  32979  32979      +       10     3     7
## 68  KK215283.1  32986  32986      +       10     3     7
## 69  KK215283.1  32989  32989      +       10     3     7
## 70  KK215283.1  32994  32994      +       10     2     8
## 71  KK215283.1  33021  33021      +       11     1    10
## 72  KK215283.1  33050  33050      +       11     2     9
## 73  KK215283.1  33081  33081      +       11     1    10
## 74  KK215283.1  33091  33091      +       12     0    12
## 75  KK215283.1  33150  33150      +       10     0    10
## 76  KK215283.1  34458  34458      +       11    10     1
## 77  KK215283.1  34484  34484      +       11     9     2
## 78  KK215283.1  34932  34932      +       10    10     0
## 79  KK215283.1  34937  34937      +       10    10     0
## 80  KK215283.1  37099  37099      +       10    10     0
## 81  KK215283.1  37831  37831      +       13    10     3
## 82  KK215283.1  38700  38700      +       10     9     1
## 83  KK215283.1  38720  38720      +       11    11     0
## 84  KK215283.1  38786  38786      +       11    11     0
## 85  KK215283.1  38856  38856      +       10    10     0
## 86  KK215283.1  39566  39566      +       10    10     0
## 87  KK215283.1  39576  39576      +       11     8     3
## 88  KK215283.1  40597  40597      +       13     4     9
## 89  KK215283.1  41681  41681      +       10    10     0
## 90  KK215283.1  42120  42120      +       10     9     1
## 91  KK215283.1  42134  42134      +       10     7     3
## 92  KK215283.1  42189  42189      +       10     9     1
## 93  KK215283.1  42294  42294      +       10     9     1
## 94  KK215283.1  42309  42309      +       14    13     1
## 95  KK215283.1  42337  42337      +       13    12     1
## 96  KK215283.1  42368  42368      +       13    11     2
## 97  KK215283.1  42371  42371      +       13    12     1
## 98  KK215283.1  42374  42374      +       11    10     1
## 99  KK215283.1  42399  42399      +       11    10     1
## 100 KK215283.1  42404  42404      +       12    10     2
## 101 KK215283.1  42407  42407      +       12    11     1
## 102 KK215283.1  42473  42473      +       12    11     1
## 103 KK215283.1  42504  42504      +       11    11     0
## 104 KK215283.1  42533  42533      +       10     7     3
## 105 KK215283.1  43073  43073      +       10    10     0
## 106 KK215283.1  43091  43091      +       11    11     0
## 107 KK215283.1  43094  43094      +       10    10     0
## 108 KK215283.1  43104  43104      +       10    10     0
## 109 KK215283.1  44555  44555      +       11     9     2
## 110 KK215283.1  44606  44606      +       11    11     0
## 111 KK215283.1  44618  44618      +       10     8     2
## 112 KK215283.1  44694  44694      +       14     0    14
## 113 KK215283.1  44710  44710      +       11     7     4
## 114 KK215283.1  44789  44789      +       10     9     1
## 115 KK215283.1  44990  44990      +       11    10     1
## 116 KK215283.1  45027  45027      +       11    10     1
## 117 KK215283.1  45074  45074      +       11    11     0
## 118 KK215283.1  45377  45377      +       10     7     3
## 119 KK215283.1  45391  45391      +       10     8     2
## 120 KK215283.1  46352  46352      +       11     5     6
## 121 KK215283.1  46403  46403      +       10     5     5
## 122 KK215283.1  46459  46459      +       10     0    10
## 123 KK215283.1  47013  47013      +       10     5     5
## 124 KK215283.1  47586  47586      +       10     7     3
## 125 KK215283.1  47672  47672      +       10     6     4
## 126 KK215283.1  47675  47675      +       10     6     4
## 127 KK215283.1  47677  47677      +       10     6     4
## 128 KK215283.1  47685  47685      +       10     6     4
## 129 KK215283.1  47690  47690      +       10     3     7
## 130 KK215283.1  47876  47876      +       10    10     0
## 131 KK215283.1  48080  48080      +       10     6     4
## 132 KK215283.1  49229  49229      +       10     6     4
## 133 KK215283.1  55651  55651      +       10     0    10
## 134 KK215283.1  56955  56955      +       10     8     2
## 135 KK215283.1  58768  58768      +       10     1     9
## 136 KK215283.1  58770  58770      +       10     0    10
## 137 KK215283.1  59050  59050      +       11     0    11
## 138 KK215283.1  59058  59058      +       10     0    10
## 139 KK215283.1  59363  59363      +       11     0    11
## 140 KK215283.1  62888  62888      +       13    13     0
## 141 KK215283.1  62927  62927      +       11    11     0
## 142 KK215283.1  63931  63931      +       12    12     0
## 143 KK215283.1  65776  65776      +       10    10     0
## 144 KK215283.1  65860  65860      +       13    11     2
## 145 KK215283.1  65896  65896      +       11     9     2
## 146 KK215283.1  65905  65905      +       11     9     2
## 147 KK215283.1  65912  65912      +       11     9     2
## 148 KK215283.1  66304  66304      +       10    10     0
## 149 KK215283.1  66963  66963      +       10     9     1
## 150 KK215283.1  67009  67009      +       12     0    12
## 151 KK215283.1  67014  67014      +       12    12     0
## 152 KK215283.1  67048  67048      +       14    14     0
## 153 KK215283.1  67118  67118      +       11    11     0
## 154 KK215283.1  67139  67139      +       10     8     2
## 155 KK215283.1  67566  67566      +       10     9     1
## 156 KK215283.1  67599  67599      +       11    10     1
## 157 KK215283.1  67605  67605      +       11    10     1
## 158 KK215283.1  67610  67610      +       11     9     2
## 159 KK215283.1  67614  67614      +       11    10     1
## 160 KK215283.1  67619  67619      +       10     8     2
## 161 KK215283.1  67623  67623      +       10     8     2
## 162 KK215283.1  67627  67627      +       10     8     2
## 163 KK215283.1  68352  68352      +       11    10     1
## 164 KK215283.1  68354  68354      +       11    10     1
## 165 KK215283.1  69173  69173      +       12    12     0
## 166 KK215283.1  69176  69176      +       12    10     2
## 167 KK215283.1  69202  69202      +       11    10     1
## 168 KK215283.1  69228  69228      +       10    10     0
## 169 KK215283.1  70298  70298      +       10     8     2
## 170 KK215283.1  70918  70918      +       12    10     2
## 171 KK215283.1  71266  71266      +       11    10     1
## 172 KK215283.1  71441  71441      +       10     9     1
## 173 KK215283.1  71472  71472      +       11     9     2
## 174 KK215283.1  71961  71961      +       11    11     0
## 175 KK215283.1  71989  71989      +       12    12     0
## 176 KK215283.1  71996  71996      +       12    12     0
## 177 KK215283.1  72001  72001      +       11    11     0
## 178 KK215283.1  72041  72041      +       10    10     0
## 179 KK215283.1  72057  72057      +       12    12     0
## 180 KK215283.1  72062  72062      +       10    10     0
## 181 KK215283.1  72863  72863      +       10     9     1
## 182 KK215283.1  74099  74099      +       11    11     0
## 183 KK215283.1  78570  78570      +       10     7     3
## 184 KK215283.1  78600  78600      +       10     6     4
## 185 KK215283.1  78650  78650      +       10     8     2
## 186 KK215283.1  78659  78659      +       12     9     3
## 187 KK215283.1  78675  78675      +       12     9     3
## 188 KK215283.1  78687  78687      +       10     7     3
## 189 KK215283.1  78823  78823      +       10     9     1
## 190 KK215283.1  78981  78981      +       10     0    10
## 191 KK215283.1  78991  78991      +       11     9     2
## 192 KK215283.1  79162  79162      +       10    10     0
## 193 KK215283.1  81657  81657      +       10    10     0
## 194 KK215283.1  81665  81665      +       10    10     0
## 195 KK215283.1  81717  81717      +       10     9     1
## 196 KK215283.1  82702  82702      +       11     9     2
## 197 KK215283.1  82727  82727      +       13    11     2
## 198 KK215283.1  82743  82743      +       12     9     3
## 199 KK215283.1  82745  82745      +       10     8     2
## 200 KK215283.1  82816  82816      +       10     8     2
## 201 KK215283.1  82929  82929      +       10     6     4
## 202 KK215283.1  82956  82956      +       12     7     5
## 203 KK215283.1  87012  87012      +       10     6     4
## 204 KK215283.1  87126  87126      +       10     8     2
## 205 KK215283.1  87238  87238      +       10     7     3
## 206 KK215283.1  87240  87240      +       10     8     2
## 207 KK215283.1  87250  87250      +       12    10     2
## 208 KK215283.1  87265  87265      +       13     9     4
## 209 KK215283.1  87325  87325      +       12    11     1
## 210 KK215283.1  87330  87330      +       12    12     0
## 211 KK215283.1  87340  87340      +       11     9     2
## 212 KK215283.1  87345  87345      +       11    11     0
## 213 KK215283.1  87356  87356      +       10    10     0
## 214 KK215283.1  87962  87962      +       11    11     0
## 215 KK215283.1  88241  88241      +       11     9     2
## 216 KK215283.1  89147  89147      +       10     7     3
## 217 KK215283.1  89152  89152      +       11     7     4
## 218 KK215283.1  89172  89172      +       11     8     3
## 219 KK215283.1  89736  89736      +       10     2     8
## 220 KK215283.1  89792  89792      +       10     2     8
## 221 KK215283.1  89800  89800      +       11     4     7
## 222 KK215283.1  89805  89805      +       11     5     6
## 223 KK215283.1  89814  89814      +       12     4     8
## 224 KK215283.1  89816  89816      +       12     6     6
## 225 KK215283.1  89965  89965      +       10     6     4
## 226 KK215283.1  90163  90163      +       11     5     6
## 227 KK215283.1  90172  90172      +       10     7     3
## 228 KK215283.1  94049  94049      +       11     7     4
## 229 KK215283.1  97557  97557      +       10     8     2
## 230 KK215283.1 104649 104649      +       10    10     0
## 231 KK215283.1 104656 104656      +       11    11     0
## 232 KK215283.1 104668 104668      +       11    10     1
## 233 KK215283.1 105036 105036      +       11    10     1
## 234 KK215283.1 105061 105061      +       11    10     1
## 235 KK215283.1 105610 105610      +       10     8     2
## 236 KK215283.1 105633 105633      +       11     9     2
## 237 KK215283.1 105657 105657      +       10    10     0
## 238 KK215283.1 108913 108913      +       11    10     1
## 239 KK215283.1 108937 108937      +       10    10     0
## 240 KK215283.1 108963 108963      +       11    11     0
## 241 KK215283.1 109818 109818      +       10     6     4
## 242 KK215283.1 110802 110802      +       10     7     3
## 243 KK215283.1 112149 112149      +       10     9     1
## 244 KK215283.1 112850 112850      +       11     7     4
## 245 KK215283.1 113246 113246      +       11    10     1
## 246 KK215283.1 113252 113252      +       11    10     1
## 247 KK215283.1 113613 113613      +       10     9     1
## 248 KK215283.1 113614 113614      +       10     7     3
## 249 KK215283.1 113619 113619      +       10     6     4
## 250 KK215283.1 114345 114345      +       12    10     2
## 251 KK215283.1 114794 114794      +       10     9     1
## 252 KK215283.1 116970 116970      +       10     9     1
## 253 KK215283.1 118884 118884      +       12    11     1
## 254 KK215283.1 118891 118891      +       10    10     0
## 255 KK215283.1 119085 119085      +       10     9     1
## 256 KK215283.1 119098 119098      +       10     9     1
## 257 KK215283.1 127585 127585      +       10    10     0
## 258 KK215283.1 131561 131561      +       10     9     1
## 259 KK215283.1 134319 134319      +       11    10     1
## 260 KK215283.1 134366 134366      +       12     9     3
## 261 KK215283.1 134373 134373      +       10     7     3
## 262 KK215283.1 135191 135191      +       13     5     8
## 263 KK215283.1 135194 135194      +       11     6     5
## 264 KK215283.1 135199 135199      +       13     8     5
## 265 KK215283.1 135217 135217      +       12     6     6
## 266 KK215283.1 135221 135221      +       11     6     5
## 267 KK215283.1 135224 135224      +       10     7     3
## 268 KK215283.1 135225 135225      +       12     8     4
## 269 KK215283.1 135256 135256      +       10     3     7
## 270 KK215283.1 135259 135259      +       10     7     3
## 271 KK215283.1 135287 135287      +       11     6     5
## 272 KK215283.1 135288 135288      +       10     2     8
## 273 KK215283.1 135294 135294      +       11     7     4
## 274 KK215283.1 135295 135295      +       11     7     4
## 275 KK215283.1 135299 135299      +       10     5     5
## 276 KK215283.1 135300 135300      +       10     5     5
## 277 KK215283.1 135305 135305      +       10     5     5
## 278 KK215283.1 135313 135313      +       10     4     6
## 279 KK215283.1 135328 135328      +       10     2     8
## 280 KK215283.1 135831 135831      +       10     4     6
## 281 KK215283.1 135914 135914      +       11     9     2
## 282 KK215283.1 137237 137237      +       10    10     0
## 283 KK215283.1 138047 138047      +       10     9     1
## 284 KK215283.1 138054 138054      +       10     9     1
## 285 KK215283.1 138059 138059      +       10     8     2
## 286 KK215283.1 138631 138631      +       10     6     4
## 287 KK215283.1 138667 138667      +       10     6     4
## 288 KK215283.1 138672 138672      +       10     3     7
## 289 KK215283.1 139876 139876      +       10    10     0
## 290 KK215283.1 139882 139882      +       10    10     0
## 291 KK215283.1 139893 139893      +       10     9     1
## 292 KK215283.1 140049 140049      +       10    10     0
## 293 KK215283.1 140051 140051      +       10    10     0
## 294 KK215283.1 140061 140061      +       11    11     0
## 295 KK215283.1 140064 140064      +       11    11     0
## 296 KK215283.1 140066 140066      +       11     7     4
## 297 KK215283.1 140076 140076      +       11    11     0
## 298 KK215283.1 140079 140079      +       11    11     0
## 299 KK215283.1 140081 140081      +       11    11     0
## 300 KK215283.1 140091 140091      +       11     7     4
## 301 KK215283.1 140096 140096      +       11     7     4
## 302 KK215283.1 140228 140228      +       10     9     1
## 303 KK215283.1 140858 140858      +       11     9     2
## 304 KK215283.1 140888 140888      +       10     9     1
## 305 KK215283.1 140902 140902      +       11    10     1
## 306 KK215283.1 140924 140924      +       14    13     1
## 307 KK215283.1 140980 140980      +       11    11     0
## 308 KK215283.1 141012 141012      +       11    11     0
## 309 KK215283.1 141086 141086      +       10    10     0
## 310 KK215283.1 141160 141160      +       11    10     1
## 311 KK215283.1 141170 141170      +       12    12     0
## 312 KK215283.1 141177 141177      +       11    10     1
## 313 KK215283.1 141184 141184      +       10    10     0
## 314 KK215283.1 141192 141192      +       10    10     0
## 315 KK215283.1 141198 141198      +       10    10     0
## 316 KK215283.1 141206 141206      +       11    11     0
## 317 KK215283.1 141211 141211      +       11    11     0
## 318 KK215283.1 142124 142124      +       10     5     5
## 319 KK215283.1 142156 142156      +       10     3     7
## 320 KK215283.1 142171 142171      +       11     5     6
## 321 KK215283.1 144017 144017      +       10     9     1
## 322 KK215283.1 144088 144088      +       10    10     0
## 323 KK215283.1 145044 145044      +       12    12     0
## 324 KK215283.1 145089 145089      +       12    10     2
## 325 KK215283.1 145104 145104      +       12    10     2
## 326 KK215283.1 147395 147395      +       10     3     7
## 327 KK215283.1 147436 147436      +       11     5     6
## 328 KK215283.1 148438 148438      +       10     9     1
## 329 KK215283.1 148452 148452      +       10    10     0
## 330 KK215283.1 148455 148455      +       10    10     0
## 331 KK215283.1 151569 151569      +       10     3     7
## 332 KK215283.1 152962 152962      +       10    10     0
## 333 KK215283.1 156476 156476      +       11     8     3
## 334 KK215283.1 156481 156481      +       11     8     3
## 335 KK215283.1 156503 156503      +       12     7     5
## 336 KK215283.1 156518 156518      +       10     7     3
## 337 KK215283.1 156882 156882      +       10     8     2
## 338 KK215283.1 156886 156886      +       10     6     4
## 339 KK215283.1 156938 156938      +       12     7     5
## 340 KK215283.1 156951 156951      +       11     9     2
## 341 KK215283.1 156967 156967      +       11     5     6
## 342 KK215283.1 156978 156978      +       12     6     6
## 343 KK215283.1 156980 156980      +       11     7     4
## 344 KK215283.1 157180 157180      +       10     8     2
## 345 KK215283.1 157183 157183      +       10     6     4
## 346 KK215283.1 157281 157281      +       10     5     5
## 347 KK215283.1 157285 157285      +       10     5     5
## 348 KK215283.1 160536 160536      +       10     7     3
## 349 KK215283.1 160576 160576      +       13     8     5
## 350 KK215283.1 160596 160596      +       11     7     4
## 351 KK215283.1 160626 160626      +       15    13     2
## 352 KK215283.1 160636 160636      +       15     6     9
## 353 KK215283.1 160646 160646      +       13    10     3
## 354 KK215283.1 162729 162729      +       10     6     4
## 355 KK215283.1 163100 163100      +       10    10     0
## 356 KK215283.1 164088 164088      +       10     9     1
## 357 KK215283.1 164100 164100      +       10     8     2
## 358 KK215283.1 164160 164160      +       10     9     1
## 359 KK215283.1 165919 165919      +       10     7     3
## 360 KK215283.1 165928 165928      +       10     5     5
## 361 KK215283.1 170217 170217      +       10    10     0
## 362 KK215283.1 170235 170235      +       11    11     0
## 363 KK215283.1 170244 170244      +       11    11     0
## 364 KK215283.1 170247 170247      +       11    11     0
## 365 KK215283.1 170255 170255      +       10    10     0
## 366 KK215283.1 170349 170349      +       12    12     0
## 367 KK215283.1 170357 170357      +       12    11     1
## 368 KK215283.1 170360 170360      +       12    12     0
## 369 KK215283.1 170366 170366      +       11    11     0
## 370 KK215283.1 170956 170956      +       10     1     9
## 371 KK215283.1 173240 173240      +       13     4     9
## 372 KK215283.1 173258 173258      +       13     7     6
## 373 KK215283.1 173270 173270      +       13     5     8
## 374 KK215283.1 173296 173296      +       13     6     7
## 375 KK215283.1 173311 173311      +       11     4     7
## 376 KK215283.1 173315 173315      +       12     5     7
## 377 KK215283.1 173327 173327      +       11     6     5
## 378 KK215283.1 173338 173338      +       11     7     4
## 379 KK215283.1 173341 173341      +       11     7     4
## 380 KK215283.1 173346 173346      +       10     6     4
## 381 KK215283.1 173349 173349      +       10     7     3
## 382 KK215283.1 173355 173355      +       11     7     4
## 383 KK215283.1 173675 173675      +       10     8     2
## 384 KK215283.1 174662 174662      +       10     9     1
## 385 KK215283.1 174749 174749      +       10    10     0
## 386 KK215283.1 177479 177479      +       10     6     4
## 387 KK215283.1 178687 178687      +       10     7     3
## 388 KK215283.1 178696 178696      +       11     7     4
## 389 KK215283.1 182396 182396      +       11     7     4
## 390 KK215283.1 182423 182423      +       10     6     4
## 391 KK215283.1 182886 182886      +       10     5     5
## 392 KK215283.1 182890 182890      +       10     7     3
## 393 KK215283.1 183275 183275      +       10     3     7
## 394 KK215283.1 187169 187169      +       10    10     0
## 395 KK215283.1 192320 192320      +       11     4     7
## 396 KK215283.1 192351 192351      +       11     6     5
## 397 KK215283.1 192396 192396      +       13     4     9
## 398 KK215283.1 192432 192432      +       12     5     7
## 399 KK215283.1 192441 192441      +       11     4     7
## 400 KK215283.1 192449 192449      +       13     7     6
## 401 KK215283.1 193447 193447      +       11     8     3
## 402 KK215283.1 193738 193738      +       10     6     4
## 403 KK215283.1 197639 197639      +       10     9     1
## 404 KK215283.1 197640 197640      +       10     9     1
## 405 KK215283.1 197644 197644      +       10    10     0
## 406 KK215283.1 197645 197645      +       11    10     1
## 407 KK215283.1 197674 197674      +       13    12     1
## 408 KK215283.1 197678 197678      +       13    13     0
## 409 KK215283.1 197683 197683      +       12    12     0
## 410 KK215283.1 197685 197685      +       12    12     0
## 411 KK215283.1 197697 197697      +       11    10     1
## 412 KK215283.1 197704 197704      +       10    10     0
## 413 KK215283.1 197720 197720      +       10    10     0
## 414 KK215283.1 199955 199955      +       10     0    10
## 415 KK215283.1 199959 199959      +       10     4     6
## 416 KK215283.1 200127 200127      +       10     1     9
## 417 KK215283.1 200180 200180      +       10     1     9
## 418 KK215283.1 200215 200215      +       10     0    10
## 419 KK215283.1 200587 200587      +       11     6     5
## 420 KK215283.1 200628 200628      +       12     7     5
## 421 KK215283.1 200639 200639      +       11     9     2
## 422 KK215283.1 200663 200663      +       11     8     3
## 423 KK215283.1 201801 201801      +       10     9     1
## 424 KK215283.1 201805 201805      +       10     9     1
## 425 KK215283.1 203495 203495      +       12    11     1
## 426 KK215283.1 203500 203500      +       11    10     1
## 427 KK215283.1 203570 203570      +       10     5     5
## 428 KK215283.1 203645 203645      +       11    11     0
## 429 KK215283.1 203724 203724      +       11    11     0
## 430 KK215283.1 203777 203777      +       10     9     1
## 431 KK215283.1 204841 204841      +       10     6     4
## 432 KK215283.1 206562 206562      +       10     1     9
## 433 KK215283.1 206637 206637      +       10     0    10
## 434 KK215283.1 206678 206678      +       11     1    10
## 435 KK215283.1 206681 206681      +       12     0    12
## 436 KK215283.1 206686 206686      +       12     0    12
## 437 KK215283.1 206703 206703      +       10     0    10
## 438 KK215283.1 206712 206712      +       11     1    10
## 439 KK215283.1 206719 206719      +       10     1     9
## 440 KK215283.1 206727 206727      +       10     0    10
## 441 KK215283.1 208819 208819      +       11    10     1
## 442 KK215283.1 208839 208839      +       10     7     3
## 443 KK215283.1 208901 208901      +       10     0    10
## 444 KK215283.1 208911 208911      +       10     5     5
## 445 KK215283.1 209476 209476      +       12     9     3
## 446 KK215283.1 209577 209577      +       10     7     3
## 447 KK215283.1 209602 209602      +       11    10     1
## 448 KK215283.1 209608 209608      +       12    11     1

Now format all the data into methylkit format

This might take a few minutes - get a cuppa.

basenames <- c("Clear2F-R1","Clear2F-R2","Clear2F-R3","Foundation","Green3F-R1","Green3F-R2","Green3F-R3",
               "Lilac4F-R1","Lilac4F-R2","Lilac4F-R3")
for (i in basenames){
  bedgz <- gsub("$","\\.bam\\.vcf\\.gz\\.cg\\.bed\\.gz",i)
  bedgz <- paste("meth_data/",bedgz,sep = "")
  mk <- gsub("$","\\.mk",i)
  mk <- paste("meth_data/",mk,sep = "")
  message(mk)
  if(file.exists(mk)){
    message(paste(mk,"exists"))
  }else{
    biscuit2methylkit(vcfz = bedgz)  
    }}
## meth_data/Clear2F-R1.mk
## meth_data/Clear2F-R1.mk exists
## meth_data/Clear2F-R2.mk
## meth_data/Clear2F-R2.mk exists
## meth_data/Clear2F-R3.mk
## meth_data/Clear2F-R3.mk exists
## meth_data/Foundation.mk
## meth_data/Foundation.mk exists
## meth_data/Green3F-R1.mk
## meth_data/Green3F-R1.mk exists
## meth_data/Green3F-R2.mk
## meth_data/Green3F-R2.mk exists
## meth_data/Green3F-R3.mk
## meth_data/Green3F-R3.mk exists
## meth_data/Lilac4F-R1.mk
## meth_data/Lilac4F-R1.mk exists
## meth_data/Lilac4F-R2.mk
## meth_data/Lilac4F-R2.mk exists
## meth_data/Lilac4F-R3.mk
## meth_data/Lilac4F-R3.mk exists
myfiles2 <- list.files("meth_data",pattern = "mk$",full.names = TRUE)
myfiles2
##  [1] "meth_data/Clear2F-R1.mk" "meth_data/Clear2F-R2.mk"
##  [3] "meth_data/Clear2F-R3.mk" "meth_data/Foundation.mk"
##  [5] "meth_data/Green3F-R1.mk" "meth_data/Green3F-R2.mk"
##  [7] "meth_data/Green3F-R3.mk" "meth_data/Lilac4F-R1.mk"
##  [9] "meth_data/Lilac4F-R2.mk" "meth_data/Lilac4F-R3.mk"
## [11] "meth_data/test-R1.mk"
AllSamples <- as.list(myfiles2[c(grep("Clear",myfiles2) , grep("Green",myfiles2), grep("Lilac",myfiles2), grep("Foundation",myfiles2))])
AllSamples
## [[1]]
## [1] "meth_data/Clear2F-R1.mk"
## 
## [[2]]
## [1] "meth_data/Clear2F-R2.mk"
## 
## [[3]]
## [1] "meth_data/Clear2F-R3.mk"
## 
## [[4]]
## [1] "meth_data/Green3F-R1.mk"
## 
## [[5]]
## [1] "meth_data/Green3F-R2.mk"
## 
## [[6]]
## [1] "meth_data/Green3F-R3.mk"
## 
## [[7]]
## [1] "meth_data/Lilac4F-R1.mk"
## 
## [[8]]
## [1] "meth_data/Lilac4F-R2.mk"
## 
## [[9]]
## [1] "meth_data/Lilac4F-R3.mk"
## 
## [[10]]
## [1] "meth_data/Foundation.mk"
myobj <-  methRead(AllSamples,
           sample.id = list("clear1","clear2","clear3","green1","green2","green3","Lilac1","Lilac2","Lilac3","Foundation"),
           assembly = "hg19",
           treatment = c(0,0,0,0,0,1,1,1,1,1),
           context = "CpG",
           mincov = 3 )
## Received list of locations.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.
## Reading file.

Simple stats

And prepare for differential methylation analysis.

lapply(myobj,function(x){
  getMethylationStats(x,plot=FALSE,both.strands=FALSE)
  getMethylationStats(x,plot=TRUE,both.strands=FALSE)
  getCoverageStats(x,plot=TRUE,both.strands=FALSE)
})
## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   60.00   83.33   73.34  100.00  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  16.66667  50.00000  66.66667  75.00000  83.33333  90.90909 100.00000 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
## 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   63.64   84.21   73.47   94.44  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  54.54545  70.00000  77.77778  84.21053  88.88889  92.30769 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
## 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   62.50   83.33   72.98   93.33  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  54.54545  68.75000  77.77778  83.33333  88.23529  91.66667 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
## 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   63.64   83.33   72.87   92.86  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  54.54545  69.23077  77.77778  83.33333  87.50000  91.66667 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
##  94.73684 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   62.50   83.33   73.35   94.12  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  54.54545  69.23077  77.77778  83.33333  88.88889  92.30769 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
## 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   62.50   82.14   72.06   92.31  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  54.54545  68.75000  76.47059  82.14286  86.66667  90.47619 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
##  93.75000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   63.64   84.21   73.39   93.33  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  55.55556  70.00000  78.57143  84.21053  88.88889  92.30769 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
##  96.15385 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   61.54   83.33   72.36   92.86  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  18.75000  53.33333  66.66667  76.92308  83.33333  87.50000  90.90909 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
## 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   62.50   83.33   73.18  100.00  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  16.66667  50.00000  66.66667  77.77778  83.33333  88.88889  94.44444 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
## 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## methylation statistics per base
## summary:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    0.00   68.18   88.00   75.45   95.45  100.00 
## percentiles:
##        0%       10%       20%       30%       40%       50%       60%       70% 
##   0.00000  20.00000  57.69231  75.00000  83.33333  88.00000  91.30435  94.44444 
##       80%       90%       95%       99%     99.5%     99.9%      100% 
##  96.15385 100.00000 100.00000 100.00000 100.00000 100.00000 100.00000

## [[1]]
## NULL
## 
## [[2]]
## NULL
## 
## [[3]]
## NULL
## 
## [[4]]
## NULL
## 
## [[5]]
## NULL
## 
## [[6]]
## NULL
## 
## [[7]]
## NULL
## 
## [[8]]
## NULL
## 
## [[9]]
## NULL
## 
## [[10]]
## NULL
myobj <- filterByCoverage(myobj,lo.count=1,lo.perc=NULL,
                                      hi.count=NULL,hi.perc=99.9)

tiles <- tileMethylCounts(myobj,win.size=1000,step.size=1000,cov.bases = 10)
meth <- unite(myobj, destrand=FALSE)
## uniting...
remove(myobj)
head(meth)
##          chr start  end strand coverage1 numCs1 numTs1 coverage2 numCs2 numTs2
## 1 KK215283.1   652  652      +         3      3      0         7      7      0
## 2 KK215283.1   772  772      +         3      3      0         9      7      2
## 3 KK215283.1   775  775      +         4      2      2         9      7      2
## 4 KK215283.1   846  846      +         4      4      0         7      6      1
## 5 KK215283.1  1662 1662      +         3      2      1        11     10      1
## 6 KK215283.1  1693 1693      +         4      2      2         7      7      0
##   coverage3 numCs3 numTs3 coverage4 numCs4 numTs4 coverage5 numCs5 numTs5
## 1        10      8      2        16     15      1         7      5      2
## 2        11     11      0        12     10      2         7      7      0
## 3        10      9      1        12      8      4         6      6      0
## 4         4      3      1         9      9      0         7      7      0
## 5         3      2      1         9      6      3         5      5      0
## 6         4      3      1        12     11      1         4      4      0
##   coverage6 numCs6 numTs6 coverage7 numCs7 numTs7 coverage8 numCs8 numTs8
## 1        16     14      2        12     12      0        13     11      2
## 2        18     17      1        15     15      0        11      9      2
## 3        18     17      1        15     14      1        12     10      2
## 4        13      9      4         8      7      1         4      3      1
## 5         8      8      0         7      6      1         6      5      1
## 6         8      7      1         9      5      4         5      4      1
##   coverage9 numCs9 numTs9 coverage10 numCs10 numTs10
## 1         8      8      0         16      15       1
## 2        11     10      1         10       9       1
## 3        10      8      2          8       6       2
## 4         3      3      0         11      11       0
## 5         4      4      0         12      11       1
## 6         5      4      1         12      12       0
# plots
getCorrelation(meth,plot=TRUE)
##               clear1    clear2    clear3    green1    green2    green3
## clear1     1.0000000 0.8010335 0.7972176 0.8089492 0.8004106 0.8102666
## clear2     0.8010335 1.0000000 0.8514886 0.8636977 0.8549608 0.8652150
## clear3     0.7972176 0.8514886 1.0000000 0.8594294 0.8506624 0.8617598
## green1     0.8089492 0.8636977 0.8594294 1.0000000 0.8631310 0.8739499
## green2     0.8004106 0.8549608 0.8506624 0.8631310 1.0000000 0.8647369
## green3     0.8102666 0.8652150 0.8617598 0.8739499 0.8647369 1.0000000
## Lilac1     0.8100390 0.8646971 0.8612553 0.8733649 0.8642651 0.8753751
## Lilac2     0.7959160 0.8496160 0.8472358 0.8580093 0.8493268 0.8606130
## Lilac3     0.7811309 0.8341844 0.8299141 0.8427774 0.8334585 0.8438452
## Foundation 0.8236716 0.8803310 0.8714570 0.8894215 0.8790812 0.8872556
##               Lilac1    Lilac2    Lilac3 Foundation
## clear1     0.8100390 0.7959160 0.7811309  0.8236716
## clear2     0.8646971 0.8496160 0.8341844  0.8803310
## clear3     0.8612553 0.8472358 0.8299141  0.8714570
## green1     0.8733649 0.8580093 0.8427774  0.8894215
## green2     0.8642651 0.8493268 0.8334585  0.8790812
## green3     0.8753751 0.8606130 0.8438452  0.8872556
## Lilac1     1.0000000 0.8602507 0.8432553  0.8865881
## Lilac2     0.8602507 1.0000000 0.8283482  0.8654853
## Lilac3     0.8432553 0.8283482 1.0000000  0.8595416
## Foundation 0.8865881 0.8654853 0.8595416  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: 10
hc = clusterSamples(meth, dist="correlation", method="ward", plot=FALSE)
## The "ward" method has been renamed to "ward.D"; note new "ward.D2"
PCASamples(meth, screeplot=TRUE)

PCASamples(meth)

Correlation heatmap

perc <- percMethylation(meth)
heatmap.2(cor(perc,method="pearson"),trace="none",scale="none",margins=c(7,7),main="Pearson")

heatmap.2(cor(perc,method="spearman"),trace="none",scale="none",margins=c(7,7),main="Spearman")

Aggregate to 1kbp tiles

#tiles <- tileMethylCounts(myobj,win.size=1000,step.size=1000,cov.bases = 10)
tiles <- unite(tiles, destrand=FALSE)
## uniting...
head(tiles)
##          chr start  end strand coverage1 numCs1 numTs1 coverage2 numCs2 numTs2
## 1 KK215283.1  1001 2000      *        91     67     24       333    302     31
## 2 KK215283.1  2001 3000      *        90     73     17       144    108     36
## 3 KK215283.1  3001 4000      *       111     97     14       336    284     52
## 4 KK215283.1  4001 5000      *       139    115     24       269    235     34
## 5 KK215283.1  5001 6000      *       463    439     24      1133   1078     55
## 6 KK215283.1  6001 7000      *       403    370     33       978    919     59
##   coverage3 numCs3 numTs3 coverage4 numCs4 numTs4 coverage5 numCs5 numTs5
## 1       230    201     29       316    278     38       214    198     16
## 2       153    126     27       180    131     49       147    107     40
## 3       338    308     30       366    326     40       348    320     28
## 4       227    204     23       308    277     31       294    257     37
## 5       812    778     34      1570   1460    110      1084    973    111
## 6      1011    904    107      1286   1184    102       756    718     38
##   coverage6 numCs6 numTs6 coverage7 numCs7 numTs7 coverage8 numCs8 numTs8
## 1       355    321     34       326    277     49       221    183     38
## 2       145    101     44       214    170     44       133     93     40
## 3       480    411     69       387    344     43       371    316     55
## 4       423    359     64       379    323     56       325    286     39
## 5      1674   1544    130      1253   1169     84      1018    913    105
## 6      1368   1240    128      1253   1162     91       699    643     56
##   coverage9 numCs9 numTs9 coverage10 numCs10 numTs10
## 1       166    151     15        559     508      51
## 2        83     61     22        314     240      74
## 3       290    255     35        677     622      55
## 4       231    198     33        640     558      82
## 5       915    847     68       2215    2139      76
## 6       865    792     73       1680    1572     108
getCorrelation(tiles,plot=TRUE)
##               clear1    clear2    clear3    green1    green2    green3
## clear1     1.0000000 0.9540509 0.9529642 0.9558521 0.9542974 0.9561720
## clear2     0.9540509 1.0000000 0.9702236 0.9733221 0.9715975 0.9734557
## clear3     0.9529642 0.9702236 1.0000000 0.9717734 0.9704185 0.9730248
## green1     0.9558521 0.9733221 0.9717734 1.0000000 0.9735249 0.9755674
## green2     0.9542974 0.9715975 0.9704185 0.9735249 1.0000000 0.9737691
## green3     0.9561720 0.9734557 0.9730248 0.9755674 0.9737691 1.0000000
## Lilac1     0.9558441 0.9731208 0.9726377 0.9750421 0.9733495 0.9758018
## Lilac2     0.9517172 0.9688151 0.9696202 0.9706265 0.9693359 0.9722097
## Lilac3     0.9495700 0.9668239 0.9649684 0.9689414 0.9668143 0.9686881
## Foundation 0.9431059 0.9612601 0.9548864 0.9636135 0.9606208 0.9597633
##               Lilac1    Lilac2    Lilac3 Foundation
## clear1     0.9558441 0.9517172 0.9495700  0.9431059
## clear2     0.9731208 0.9688151 0.9668239  0.9612601
## clear3     0.9726377 0.9696202 0.9649684  0.9548864
## green1     0.9750421 0.9706265 0.9689414  0.9636135
## green2     0.9733495 0.9693359 0.9668143  0.9606208
## green3     0.9758018 0.9722097 0.9686881  0.9597633
## Lilac1     1.0000000 0.9720180 0.9681887  0.9591494
## Lilac2     0.9720180 1.0000000 0.9635826  0.9485543
## Lilac3     0.9681887 0.9635826 1.0000000  0.9578280
## Foundation 0.9591494 0.9485543 0.9578280  1.0000000

clusterSamples(tiles, 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: 10
hc = clusterSamples(tiles, dist="correlation", method="ward", plot=FALSE)
## The "ward" method has been renamed to "ward.D"; note new "ward.D2"
PCASamples(tiles, screeplot=TRUE)

PCASamples(tiles)

perc <- percMethylation(tiles)
heatmap.2(cor(perc,method="pearson"),trace="none",scale="none",margins=c(7,7),main="Pearson")

heatmap.2(cor(perc,method="spearman"),trace="none",scale="none",margins=c(7,7),main="Spearman")

Session information

For reproducibility

sessionInfo()
## R version 4.0.1 (2020-06-06)
## Platform: x86_64-conda_cos6-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.2 LTS
## 
## Matrix products: default
## BLAS/LAPACK: /ceph-g/opt/miniconda3/envs/R40/lib/libopenblasp-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] stats4    parallel  stats     graphics  grDevices utils     datasets 
## [8] methods   base     
## 
## other attached packages:
##  [1] seqinr_4.2-8         gplots_3.1.1         methylKit_1.14.2    
##  [4] limma_3.44.3         GenomicRanges_1.40.0 GenomeInfoDb_1.24.2 
##  [7] IRanges_2.22.2       S4Vectors_0.26.1     BiocGenerics_0.34.0 
## [10] RColorBrewer_1.1-2   dplyr_1.0.0          reshape2_1.4.4      
## [13] R.utils_2.11.0       R.oo_1.24.0          R.methodsS3_1.8.1   
## 
## loaded via a namespace (and not attached):
##  [1] mclust_5.4.7                Rcpp_1.0.7                 
##  [3] bdsmatrix_1.3-4             mvtnorm_1.1-2              
##  [5] lattice_0.20-41             Rsamtools_2.4.0            
##  [7] Biostrings_2.56.0           gtools_3.9.2               
##  [9] digest_0.6.25               R6_2.4.1                   
## [11] plyr_1.8.6                  emdbook_1.3.12             
## [13] evaluate_0.14               coda_0.19-4                
## [15] ggplot2_3.3.5               pillar_1.4.4               
## [17] zlibbioc_1.34.0             rlang_0.4.11               
## [19] data.table_1.12.8           Matrix_1.2-18              
## [21] bbmle_1.0.24                qvalue_2.20.0              
## [23] rmarkdown_2.3               splines_4.0.1              
## [25] BiocParallel_1.22.0         stringr_1.4.0              
## [27] fastseg_1.34.0              RCurl_1.98-1.5             
## [29] munsell_0.5.0               DelayedArray_0.14.1        
## [31] compiler_4.0.1              numDeriv_2016.8-1.1        
## [33] rtracklayer_1.48.0          xfun_0.15                  
## [35] pkgconfig_2.0.3             htmltools_0.5.0            
## [37] SummarizedExperiment_1.18.2 tidyselect_1.1.0           
## [39] tibble_3.0.1                GenomeInfoDbData_1.2.3     
## [41] matrixStats_0.61.0          XML_3.99-0.8               
## [43] crayon_1.3.4                GenomicAlignments_1.24.0   
## [45] MASS_7.3-51.6               bitops_1.0-7               
## [47] grid_4.0.1                  gtable_0.3.0               
## [49] lifecycle_0.2.0             magrittr_1.5               
## [51] scales_1.1.1                KernSmooth_2.23-17         
## [53] stringi_1.4.6               XVector_0.28.0             
## [55] ellipsis_0.3.1              generics_0.0.2             
## [57] vctrs_0.3.1                 tools_4.0.1                
## [59] ade4_1.7-18                 Biobase_2.48.0             
## [61] glue_1.4.1                  purrr_0.3.4                
## [63] yaml_2.2.1                  colorspace_1.4-1           
## [65] caTools_1.18.2              knitr_1.28