[BioC] DiffBind error

Rory Stark Rory.Stark at cruk.cam.ac.uk
Wed Jan 30 16:09:19 CET 2013


Hi again-

You need to explicitly retrieve the results of the analysis
(differentially bound sites) using dba.report:

> chippy.db = dba.report(chippy)

By default this will come back as a GRanges object. You can get it back as
a data frame using the DataType parameter:

> chippy.db = dba.report(chippy, DataType = DBA_DATA_FRAME)


(note you can set the default data type by setting
chippy$config$DataType=DBA_DATA_FRAME).

You can save the report to a csv file directly by specifying a file name
to dba.report:

> chippy.db = dba.report(chippy,file="ChIPseq_diffbind")

Which will return the report as a Granges object and create a file called
"ChIPseq_diffbind.csv".

Note that you can control the contents of the report with other parameters
to dba.report. By default, it will include sites with an FDR < 0.10; for
each site you will get the interval, the log2 mean concentration of
normalized reads overall and in the two groups, the log2 fold change
(positive or negative depending on which group has a greater read
concentration), as well as the p-value and FDR (corrected for multiple
testing). You can change the threshold using the "th" and "fold"
parameters (setting th=1 will include all sites in the report, sorted by
lowest FDR). You can also include the counts for each sample (normalized
or non-normalized) by setting bCounts=TRUE, among other options.

I usually return all the sites as a csv file, with counts for each sample,
to the biologist and let them think about cutoffs, using:

> chippy.db = dba.report(chippy, th=1, bCounts=T, file"dba_chippy")


Cheers-
Rory


On 30/01/2013 14:51, "jluis.lavin at unavarra.es" <jluis.lavin at unavarra.es>
wrote:

>Hello again,
>
>Thank you very much to Rory and Gordon for your kind and accurate help.
>Changing the minMembers parameter everything seemed to work fine and I've
>been able to perform the next steps of the whole analysis.
>
>Now I'm struggling with the retrieval of the analysis data; when I try to
>write the results into a .csv table with a simple command I get the
>following error:
>
>> write.csv( chippy, file="ChIPseq_diffbind.csv" )
>
>
>Error in as.data.frame.default(x[[i]], optional = TRUE, stringsAsFactors =
>stringsAsFactors) :
>  cannot coerce class '"DBA"' into a data.frame
>
>I read DiffBind's vignette and manual and I found the "save" command, but
>it doesn't allow me to retrieve the data either.
>Is there a way to coerce a DBA class object into a dataframe implemented
>in DiffBind?
>The person I'll be analyzing ChIPseq experiments for will need tables with
>differentially bound sites from his data comparisons so I need to learn
>how to get data out of DiffBind...
>
>Thanks in advance
>
>JL
>
>> sessionInfo()
>R version 2.15.1 (2012-06-22)
>Platform: x86_64-redhat-linux-gnu (64-bit)
>
>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=C                 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  stats     graphics  grDevices utils     datasets  methods
>[8] base
>
>other attached packages:
>[1] DESeq_1.10.1         lattice_0.20-10      locfit_1.5-8
>[4] DiffBind_1.4.1       Biobase_2.18.0       GenomicRanges_1.10.5
>[7] IRanges_1.16.4       BiocGenerics_0.4.0   BiocInstaller_1.8.3
>
>loaded via a namespace (and not attached):
> [1] amap_0.8-7           annotate_1.36.0      AnnotationDbi_1.20.3
> [4] DBI_0.2-5            edgeR_3.0.8          gdata_2.12.0
> [7] genefilter_1.40.0    geneplotter_1.36.0   gplots_2.11.0
>[10] grid_2.15.1          gtools_2.7.0         KernSmooth_2.23-8
>[13] limma_3.14.4         RColorBrewer_1.0-5   RSQLite_0.11.2
>[16] splines_2.15.1       stats4_2.15.1        survival_2.36-14
>[19] tools_2.15.1         XML_3.95-0.1         xtable_1.7-0
>[22] zlibbioc_1.4.0
>
>
>El Mie, 30 de Enero de 2013, 13:28, Rory Stark escribió:
>>
>> Hello-
>>
>> It looks like your sample sheet is fine.
>>
>> By default, if no contrasts are set up using dba.contrast,  when you
>>call
>> dba.analyze it attempts to add all the two-way contrasts between groups
>> where there are at least three samples in each group. In your case,
>> nothing meets these conditions, as there are only two "Resistant"
>>samples.
>>
>> If you add a call to dba.contrast before the call to dba.analyze, you
>>will
>> get the contrast you desire:
>>
>>> chippy = dba.contrast(chippy, minMembers=2)
>>
>> or, slightly more explicitly:
>>
>>> chippy = dba.contrast(chippy, categories=DBA_CONDITION, minMembers=2)
>>
>> or, even more explicitly:
>>
>>> chippy = dba.contrast(chippy, group1=chippy$masks$Resistant, group2 =
>>>chippy$masksResponsive, name1="Resistant", name2="Responsive")
>>
>> then:
>>
>>> chippy = dba.analyze(chippy)
>>
>> Cheers-
>> Rory
>>
>Hi,
>
>By default, dba.contrast won't create contrasts with less than 3 members
>in each group.  The easiest solution is to set minMembers=2 when calling
>dba.contrast:
>
>> chippy = dba.contrast(chippy, categories=DBA_CONDITION, minMembers=2)
>
>Or you can create the contrast explicitly:
>
>> chippy = dba.contrast(chippy, group1=chippy$masks$Resistant,
>>group2=chippy$masks$Responsive, name1="a name", name2="another name")
>
>Hope this helps...
>
>Cheers,
>
> - Gord
>>
>> On Tue, 29 Jan 2013 15:21:23 +0100, jluis.lavin at unavarra.es wrote:
>>
>>
>>>
>>>Dear list,
>>>
>>>I made a new samplesheet to use with DiffBind with this format:
>>>
>>>SampleID	Tissue	Factor	Condition	Replicate	bamReads	bamControl	Peaks
>>>Contrl1	Neural	K9	Resistant	1	Control1.bed	Contro_input.bed	Control1_pea
>>>ks
>>>.bed
>>>Contrl2	Neural	K9	Resistant	2	Control2.bed	Control_input.bed	Control2_pe
>>>ak
>>>s.bed
>>>A4_1	Neural	K9	Responsive	1	A4_1.bed	A4_input.bed	A4_1_peaks.bed
>>>A4_2	Neural	K9	Responsive	2	A4_2.bed	A4_input.bed	A4_2_peaks.bed
>>>A21_1	Neural	K9	Responsive	1	A21_1.bed	A21_input.bed	A21_1_peaks.bed
>>>A21_2	Neural	K9	Responsive	2	A21_2.bed	A21_input.bed	A21_2_peaks.bed
>>>
>>>I can load the files,
>>>
>>>> chippy = dba(sampleSheet="Peaksets_sample_sheet.csv")
>>>
>>>plot them
>>>
>>>>plot(chippy)
>>>
>>>and count the reads,
>>>
>>>>chippy = dba.count(chippy, minOverlap=3)
>>>
>>>but when I try to establish a contrast based on the condition metadata,
>>>I
>>>get the following warning message:
>>>
>>>> chippy = dba.contrast(chippy, categories=DBA_CONDITION)
>>>
>>>Warning message:
>>>No contrasts added. Perhaps try more categories, or lower value for
>>>minMembers.
>>>
>>>So when I try to perform the analysis, it doesn't work:
>>>
>>>> chippy = dba.analyze(chippy)
>>>
>>>Error in pv.DBA(DBA, method, bSubControl, bFullLibrarySize, bTagwise =
>>>bTagwise,  :
>>>  Unable to perform analysis: no contrasts specified.
>>>In addition: Warning message:
>>>No contrasts added. Perhaps try more categories, or lower value for
>>>minMembers.
>>>
>>>->My questions are:
>>>-What is the contrast for DiffBind (I added the input for each set of
>>>samples, and 2 biological replicates as control)?
>>>-Is there something wrong with the sample sheet?
>>>-Shouldn't the files to analyze be in bed format?
>>>
>>>Thanks in advance
>>>
>>>JL
>>>
>>>We may monitor all incoming and outgoing emails in line with current
>>>legislation. We have taken steps to ensure that this email and
>>>attachments are free from any virus, but it remains your responsibility
>>>to ensure that viruses do not adversely affect you.
>>>Cancer Research UKRegistered charity in England and Wales (1089464),
>>>Scotland (SC041666) and the Isle of Man (1103)
>>>A company limited by guarantee.  Registered company in England and Wales
>>>(4325234) and the Isle of Man (5713F).
>>>Registered Office Address: Angel Building, 407 St John Street, London
>>>EC1V 4AD.
>>
>>
>
>
>-- 
>Dr. José Luis Lavín Trueba
>
>Dpto. de Producción Agraria
>Grupo de Genética y Microbiología
>Universidad Pública de Navarra
>31006 Pamplona
>Navarra
>SPAIN
>
>



More information about the Bioconductor mailing list