[BioC] HTqPCR
Heidi Dvinge
Heidi.Dvinge at cancer.org.uk
Wed Jun 20 10:45:47 CEST 2012
Hi SIlvia,
On 19 Jun 2012, at 18:48, Silvia Halim wrote:
> Hi Heidi,
>
> Thanks for your tips. I figured I could probably use plotCtVariation(). I am able to use the function to plot variation across samples but how can I use it to plot across genes (or features)?
> I tried following commands:
> plotCtVariation(temp[,1:10], variation = "sd", log = TRUE, main = "SD of replicated features", col = "lightgrey")
> plotCtVariation(temp[1:10,], variation = "sd", log = TRUE, main = "SD of replicated features", col = "lightgrey")
Here you're just subsetting your qPCRset object before plotting, but you're not changing the actual plots.
> There's a difference in the plots but both plots give me same labels on x-axis, i.e. sample names, though I was expecting the second command would give me gene names on x-axis label.
>
If you look at the plotCtVariation help files (especially the 'Examples' and 'Details' section), the parameters sample.reps and feature.reps controls whether you plot the variation for each gene across or within samples. In order to get gene names, you have to set sample.reps, to indicate which samples are replicates of each other.
Per default, the function calculates the variation between replicated features within each of your samples, and plots the distribution (boxplot) of this variation for each sample. If you want to check individual features or samples more specifically, you ahve to use type="detail" and possibly add.featurenames=TRUE. There are some examples included in the plotCtVariation help file.
> Also, the manual says we can exclude unreliable or undetermined data by setting the Ct values to NA using filterCategory. I am wondering how I can get rid of NA data from the plate. I also cannot exclude this kind of data or those having 'Failed' flags the very first time before reading in the input as a qPCRset object because the input has to be something like 48 x 48 or 96 x 96.
>
The question is, why do you want to remove the NA values? If you just leave them as NAs, then they're ignored during e.g. the calculation of differential expression and for most plotting purposes.
You can't remove them as such, since (as you note), the object has to be in a certain features x samples format. If you want, you can replace them though, if you e.g. want to set all NA values to Ct=40: exprs(temp)[is.na(exprs(temp))] <- 40. But beware, because in that case the value '40' will be include into all numerical calculations, which may not be what you want.
HTH
\Heidi
> Many thanks,
> Silvia
>
> -----Original Message-----
> From: Heidi Dvinge
> Sent: 19 June 2012 2:05 PM
> To: Silvia Halim
> Cc: bioconductor at r-project.org
> Subject: Re: HTqPCR
>
> Hi Silvia,
>
> On 18 Jun 2012, at 17:51, Silvia Halim wrote:
>
>> Hi Heidi,
>>
>> The function breaks at plotCtReps.
>>> traceback()
>> 1: plotCtReps(temp, card = 2, percent = 20, xlim = c(0, 100), ylim = c(0,
>> 100))
>>
>> You've pointed out the problem about the duplicates as I have 3 replicates on my assay. I got confused reading the manual as it says plotCtReps can be used for a sample containing duplicate measurements (which I thought to be 2 or more measurements).
>>> table(table(featureNames(temp)))
>>
>> 3 6
>> 30 1
>>
> If you try running the examples for plotCtReps, you'll see that the function directly plots two replicates of a feature against each other on the (x,y) axis. 3D (x,y,z) plots aren't implemented, so features that are replicated 3 times can't be plotted. I'll try to clarify the text for the function.
>
> Perhaps something like plotCtVariation() will give you what you're after? If you only want to visually inspect your data, then grep("plot", ls("package:HTqPCR"), value=TRUE) will list all the plotting functions available in HTqPCR.
>
> HTH
> \Heidi
>
>> Btw there's no NA in my data.
>>> sum(is.na(temp))
>> [1] 0
>> Warning message:
>> In is.na(temp) : is.na() applied to non-(list or vector) of type 'S4'
>>>
>>
>> Thanks,
>> Silvia
>>
>> -----Original Message-----
>> From: Heidi Dvinge
>> Sent: 15 June 2012 9:06 PM
>> To: Silvia Halim
>> Cc: bioconductor at r-project.org
>> Subject: Re: HTqPCR
>>
>> Hi Silvia,
>>
>> On 15 Jun 2012, at 18:45, Silvia Halim wrote:
>>
>>> Hi Heidi,
>>>
>>> I ran into below problem when using plotCtReps.
>>>
>>>> plotCtReps(temp, card = 1, percent = 20, xlim = c(0,50), ylim =
>>>> c(0,50))
>>> Error in split.data[[s]] : subscript out of bounds In addition:
>>> Warning messages:
>>> 1: In min(x, na.rm = na.rm) :
>>> no non-missing arguments to min; returning Inf
>>> 2: In max(x, na.rm = na.rm) :
>>> no non-missing arguments to max; returning -Inf
>>>> plotCtReps(temp, card = 1, percent = 20, xlim = c(0,50), ylim =
>>>> c(0,50))
>>> Error in split.data[[s]] : subscript out of bounds In addition:
>>> Warning messages:
>>> 1: In min(x, na.rm = na.rm) :
>>> no non-missing arguments to min; returning Inf
>>> 2: In max(x, na.rm = na.rm) :
>>> no non-missing arguments to max; returning -Inf
>>>> plotCtReps(temp, card = 2, percent = 20, xlim = c(0,100), ylim =
>>>> c(0,100))
>>> Error in split.data[[s]] : subscript out of bounds In addition:
>>> Warning messages:
>>> 1: In min(x, na.rm = na.rm) :
>>> no non-missing arguments to min; returning Inf
>>> 2: In max(x, na.rm = na.rm) :
>>> no non-missing arguments to max; returning -Inf
>>
>> What's the output from traceback(), i.e. exactly where does the function break?
>>>
>> A couple of things you can try:
>>
>> - plotCtReps is meant to be used in cases where there are exactly 2 replicates of the features on your assay. Is this the case? For example, with the data below there are 190 features that will be plotted, and 1 that will be skipped:
>>> data(qPCRraw)
>>> table(table(featureNames(qPCRraw)))
>> 2 4
>> 190 1
>>
>> - are there any NAs in your data? E.g. sum(is.na(qPCRraw))>0.
>>
>> HTH
>> \Heidi
>>
>>> Here is how 'temp' looks like
>>>> temp
>>> An object of class "qPCRset"
>>> Size: 96 features, 96 samples
>>> Feature types: Reference, Test
>>> Feature names: b-Actin b-Actin b-Actin ...
>>> Feature classes:
>>> Feature categories: OK
>>> Sample names: NTC_4 PMPT352 NTC_3 ...
>>>
>>> Do you know why it is complaining about split.data?
>>>
>>> Thanks,
>>> Silvia
>>>
>>> -----Original Message-----
>>> From: Heidi Dvinge
>>> Sent: 11 June 2012 6:11 PM
>>> To: Silvia Halim
>>> Subject: Re: HTqPCR
>>>
>>> Ok, so you already have a 96 by 96 matrix, so you don't need changeCtLayout.
>>> Good luck with the rest, and let me know if you encounter any problems.
>>>
>>> On 11 Jun 2012, at 19:05, Silvia Halim wrote:
>>>
>>>> Hi Heidi,
>>>>
>>>> Thank you for your clarification.
>>>>
>>>> Btw this is how it looks like when I type 'temp'
>>>>> temp
>>>> An object of class "qPCRset"
>>>> Size: 96 features, 96 samples
>>>> Feature types: Reference, Test
>>>> Feature names: b-Actin b-Actin b-Actin ...
>>>> Feature classes:
>>>> Feature categories: OK
>>>> Sample names: NTC_4 PMPT352 NTC_3 ...
>>>>
>>>> Cheers,
>>>> Silvia
>>>>
>>>> -----Original Message-----
>>>> From: Heidi Dvinge
>>>> Sent: 08 June 2012 7:12 PM
>>>> To: Silvia Halim
>>>> Subject: Re: HTqPCR
>>>>
>>>> Hi Silvia,
>>>>
>>>> what are the dimensions of the "temp" object that you read in? I.e.
>>>> what does it look like if you just type
>>>>> temp
>>>>
>>>> If you read in the data with n.features=96 and n.data=96, then you should already have an object with 96 rows and 96 columns, in which case you don't need to change the layout.
>>>>
>>>> Best,
>>>> \Heidi
>>>>
>>>> On 8 Jun 2012, at 19:13, Silvia Halim wrote:
>>>>
>>>>> Hi Heidi,
>>>>>
>>>>> I finally have time to try out your HTqPCR bioconductor package again and I was trying to use 'changeCtLayout' function. However, I got following error message:
>>>>>
>>>>>> qPCRnew <- changeCtLayout(temp, sample.order = sample_order)
>>>>> Error in data.frame(..., check.names = FALSE) :
>>>>> arguments imply differing number of rows: 0, 96 In addition:
>>>>> Warning
>>>>> message:
>>>>> In split.default(x = seq_len(nrow(x)), f = f, drop = drop, ...) :
>>>>> data length is not a multiple of split variable
>>>>>
>>>>> The commands that I run are as follows:
>>>>>> temp <- readCtData("110614 BENIGN_1 DATA 96X96.csv", path =
>>>>>> getwd(), n.features = 96, n.data=96, flag = 9, feature = 5, type=
>>>>>> 6, Ct = 7, position = 1, skip = 12, sep = ",") sample_order <-
>>>>>> rep(sampleNames(temp), each = 96) qPCRnew <- changeCtLayout(temp,
>>>>>> sample.order = sample_order)
>>>>>
>>>>> I've tried to follow what's written in changeCtLayout function description. Can you please advise what went wrong?
>>>>>
>>>>> Thanks,
>>>>> Silvia
>>>>>
>>>>> -----Original Message-----
>>>>> From: Heidi Dvinge
>>>>> Sent: 29 April 2012 8:18 PM
>>>>> To: Silvia Halim
>>>>> Subject: Re: HTqPCR
>>>>>
>>>>> HI Silvia,
>>>>>
>>>>> I'm glad you got it working. Depending on what you're supposed to do with the data, you may need to tweak some functions slightly, as you mention. Let me know if you run into any more trouble.
>>>>>
>>>>> Cheers
>>>>> \Heidi
>>>>>
>>>>> On 26 Apr 2012, at 18:37, Silvia Halim wrote:
>>>>>
>>>>>> Hi Heidi,
>>>>>>
>>>>>> Thanks for the help! It's working for me now. Right now I'm figuring it out how I can use the functions that you described in the vignette. I might have to tweak the parameters for using the functions on Fluidigm data.
>>>>>>
>>>>>> Cheers,
>>>>>> Silvia
>>>>>>
>>>>>> -----Original Message-----
>>>>>> From: Heidi Dvinge
>>>>>> Sent: 25 April 2012 8:56 AM
>>>>>> To: Silvia Halim
>>>>>> Subject: Re: HTqPCR
>>>>>>
>>>>>> Hiya,
>>>>>>
>>>>>> sorry, I only just now realised that you'd attached a file. When I saved as csv, the following command worked:
>>>>>>
>>>>>>> raw <- readCtData("110614 BENIGN_1 DATA 96x96.csv",
>>>>>>> format="BioMark",
>>>>>>> n.features=96*96) raw
>>>>>> An object of class "qPCRset"
>>>>>> Size: 9216 features, 1 samples
>>>>>> Feature types:
>>>>>> Feature names: b-Actin b-Actin b-Actin ...
>>>>>> Feature classes:
>>>>>> Feature categories: OK
>>>>>> Sample names: 110614 BENIGN_1 DATA 96x96 ...
>>>>>>
>>>>>> The data isn't transformed into a 96x96 format immediately though (in case you read in multiple arrays, and want to normalise them independently). If you want to change this, you can use changeCtLayout(). Alternatively you can say:
>>>>>>
>>>>>>> raw <- readCtData("110614 BENIGN_1 DATA 96x96.csv",
>>>>>>> format="BioMark", n.features=96, n.data=96) raw
>>>>>> An object of class "qPCRset"
>>>>>> Size: 96 features, 96 samples
>>>>>> Feature types:
>>>>>> Feature names: b-Actin b-Actin b-Actin ...
>>>>>> Feature classes:
>>>>>> Feature categories: OK
>>>>>> Sample names: Sample1 Sample2 Sample3 ...
>>>>>>> plotCtArray(raw)
>>>>>>
>>>>>> HTH
>>>>>> \Heidi
>>>>>>
>>>>>> On 24 Apr 2012, at 17:55, Silvia Halim wrote:
>>>>>>
>>>>>>> Hi Heidi,
>>>>>>>
>>>>>>> I have some problems updating R on lustre. Therefore, I chose to run HTqPCR on my desktop for the moment.
>>>>>>>
>>>>>>> Reading in your sample file looks fine, however, reading in the
>>>>>>> file that I showed you just now gave me below error message. (The
>>>>>>> file is as attached)
>>>>>>>
>>>>>>>> temp <- readCtData("110614 BENIGN_1 DATA 96x96.xlsx", path =
>>>>>>>> getwd() , n.features = 96*96, flag = 9, feature = 5, type= 6, Ct
>>>>>>>> = 7,position = 1, skip = 12, sep = ",")
>>>>>>> Error in read.table(file = file, header = header, sep = sep, quote = quote, :
>>>>>>> no lines available in input
>>>>>>> In addition: Warning message:
>>>>>>> In readLines(file, skip) :
>>>>>>> incomplete final line found on 'C:/Users/halim01/Documents/20110627_RossAdamsH_DN_Fluid/110614 BENIGN_1 DATA 96x96.xlsx'
>>>>>>>> sessionInfo()
>>>>>>> R version 2.14.0 (2011-10-31)
>>>>>>> Platform: x86_64-pc-mingw32/x64 (64-bit)
>>>>>>>
>>>>>>> locale:
>>>>>>> [1] LC_COLLATE=English_United Kingdom.1252 LC_CTYPE=English_United Kingdom.1252 LC_MONETARY=English_United Kingdom.1252 LC_NUMERIC=C LC_TIME=English_United Kingdom.1252
>>>>>>>
>>>>>>> attached base packages:
>>>>>>> [1] stats graphics grDevices utils datasets methods base
>>>>>>>
>>>>>>> other attached packages:
>>>>>>> [1] Biostrings_2.22.0 IRanges_1.12.6 BiocInstaller_1.2.1 marray_1.32.0 HTqPCR_1.8.0 limma_3.10.3 RColorBrewer_1.0-5 Biobase_2.14.0 gdata_2.8.2
>>>>>>>
>>>>>>> loaded via a namespace (and not attached):
>>>>>>> [1] affy_1.32.1 affyio_1.22.0 gplots_2.10.1 gtools_2.6.2 preprocessCore_1.16.0 tools_2.14.0 zlibbioc_1.0.1
>>>>>>>>
>>>>>>>
>>>>>>> I did a quick check on the file and it only has 9228 lines including 12 header lines which I had skipped when reading in the file. Do you know what could possibly go wrong?
>>>>>>>
>>>>>>> Cheers,
>>>>>>> Silvia
>>>>>>>
>>>>>>> -----Original Message-----
>>>>>>> From: Heidi Dvinge
>>>>>>> Sent: 24 April 2012 5:09 PM
>>>>>>> To: Silvia Halim
>>>>>>> Subject: Re: HTqPCR
>>>>>>>
>>>>>>> Hm, that looks like it may be x11 acting up. I often have similar issues when I work on a remote server.
>>>>>>>
>>>>>>> Actually, the processing of Fluidigm files is very computationally light. So you can easily do it on your desktop, if you can't update on lustre.
>>>>>>>
>>>>>>> I can also email you and older version of the vignette if you want to have a look. However, in HTqPCR 1.2.0 I don't even think I had a dedicated function for plotting the Fluidigm assays yet (the plotCtArray shown in the vignette).
>>>>>>>
>>>>>>> Cheers
>>>>>>> \Heidi
>>>>>>>
>>>>>>> On 24 Apr 2012, at 16:39, Silvia Halim wrote:
>>>>>>>
>>>>>>>> Hi Heidi,
>>>>>>>>
>>>>>>>> This is what I got when accessing the vignette.
>>>>>>>>
>>>>>>>>> openVignette(package="HTqPCR")
>>>>>>>> Please select a vignette:
>>>>>>>>
>>>>>>>> 1: HTqPCR - qPCR analysis in R
>>>>>>>>
>>>>>>>> Selection: 1
>>>>>>>> Opening
>>>>>>>> /home/mib-cri/local/lib64/R/library/HTqPCR/doc/HTqPCR.pdf
>>>>>>>>> xprop: unable to open display ''
>>>>>>>> /usr/local/bin/xdg-open: line 370: firefox: command not found
>>>>>>>> /usr/local/bin/xdg-open: line 370: mozilla: command not found
>>>>>>>> /usr/local/bin/xdg-open: line 370: netscape: command not found
>>>>>>>> xdg-open: no method available for opening '/home/mib-cri/local/lib64/R/library/HTqPCR/doc/HTqPCR.pdf'
>>>>>>>>
>>>>>>>> Sorry for the confusion, you are right that I was looking at a newer version of HTqPCR than the one installed on lustre. I think that's because I have different installations of HTqPCR on lustre and on my desktop. If I can update the one on lustre, I'll go ahead with the update.
>>>>>>>>
>>>>>>>> Thank you,
>>>>>>>> Silvia
>>>>>>>>
>>>>>>>> -----Original Message-----
>>>>>>>> From: Heidi Dvinge
>>>>>>>> Sent: 24 April 2012 4:28 PM
>>>>>>>> To: Silvia Halim
>>>>>>>> Subject: Re: HTqPCR
>>>>>>>>
>>>>>>>> Ah, right, it looks like you have an older version of R, and therefore also HTqPCR.
>>>>>>>>
>>>>>>>> The most current release version is 1.10.0. In that version, readCtData() was modified to accept different types of input data, including from Fluidigm. Before that, this sort of data had to be read in 'manually'.
>>>>>>>>
>>>>>>>> I guess the vignette that you were looking at comes from a
>>>>>>>> version of HTqPCR that's newer than the one you have installed?
>>>>>>>> If you access the vignette corresponding to your HTqPCR version
>>>>>>>> via
>>>>>>>>> openVignette(package="HTqPCR")
>>>>>>>> what do you get then?
>>>>>>>>
>>>>>>>> If you get an older version, then depending on how old it is, there may be a section towards the end giving an example of how to process Fluidigm data more 'manually'. If not, an update may be your best bet.
>>>>>>>>
>>>>>>>> Cheers
>>>>>>>> \Heidi
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> On 24 Apr 2012, at 16:15, Silvia Halim wrote:
>>>>>>>>
>>>>>>>>> Hi Heidi,
>>>>>>>>>
>>>>>>>>> Thanks for looking into the matter. Below is the output of my
>>>>>>>>> sessionInfo()
>>>>>>>>>
>>>>>>>>>> sessionInfo()
>>>>>>>>> R version 2.13.0 (2011-04-13)
>>>>>>>>> Platform: x86_64-unknown-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=C LC_MESSAGES=en_US.UTF-8
>>>>>>>>> [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
>>>>>>>>> [9] LC_ADDRESS=C LC_TELEPHONE=C
>>>>>>>>> [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
>>>>>>>>>
>>>>>>>>> attached base packages:
>>>>>>>>> [1] stats graphics grDevices utils datasets methods base
>>>>>>>>>
>>>>>>>>> other attached packages:
>>>>>>>>> [1] marray_1.26.0 Biostrings_2.20.1 IRanges_1.10.3 HTqPCR_1.2.0
>>>>>>>>> [5] limma_3.6.9 RColorBrewer_1.0-2 Biobase_2.12.1 gdata_2.8.0
>>>>>>>>>
>>>>>>>>> loaded via a namespace (and not attached):
>>>>>>>>> [1] affy_1.26.1 affyio_1.20.0 gplots_2.8.0
>>>>>>>>> [4] gtools_2.6.2 preprocessCore_1.14.0
>>>>>>>>>>
>>>>>>>>>
>>>>>>>>> Cheers,
>>>>>>>>> Silvia
>>>>>>>>>
>>>>>>>>> -----Original Message-----
>>>>>>>>> From: Heidi Dvinge
>>>>>>>>> Sent: 24 April 2012 4:07 PM
>>>>>>>>> To: Silvia Halim
>>>>>>>>> Subject: HTqPCR
>>>>>>>>>
>>>>>>>>> Hi Silvia,
>>>>>>>>>
>>>>>>>>> I just tested the read fluidigm from the vignette, and it works on both my mac and a single unix system that I've tested. Although from the errors you were getting, it seemed like the headers weren't been read correctly/at all.
>>>>>>>>>
>>>>>>>>> Would you mind sending me the output of your sessionInfo(), so I can compare which package versions we have?
>>>>>>>>>
>>>>>>>>> Best,
>>>>>>>>> \Heidi
>>>>>>>>>
>>>>>>>>>> sessionInfo()
>>>>>>>>> R version 2.15.0 (2012-03-30)
>>>>>>>>> Platform: x86_64-apple-darwin9.8.0/x86_64 (64-bit)
>>>>>>>>>
>>>>>>>>> locale:
>>>>>>>>> [1]
>>>>>>>>> en_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8
>>>>>>>>>
>>>>>>>>> attached base packages:
>>>>>>>>> [1] tools stats graphics grDevices utils datasets methods base
>>>>>>>>>
>>>>>>>>> other attached packages:
>>>>>>>>> [1] HTqPCR_1.10.0 limma_3.12.0 RColorBrewer_1.0-5 Biobase_2.16.0
>>>>>>>>> [5] BiocGenerics_0.2.0
>>>>>>>>>
>>>>>>>>> loaded via a namespace (and not attached):
>>>>>>>>> [1] affy_1.34.0 affyio_1.24.0 BiocInstaller_1.4.3
>>>>>>>>> [4] gdata_2.8.2 gplots_2.10.1 gtools_2.6.2
>>>>>>>>> [7] preprocessCore_1.18.0 stats4_2.15.0 zlibbioc_1.2.0
>>>>>>>>>
>>>>>>>>>
>>>>>>>>
>>>>>>>
>>>>>>> <110614 BENIGN_1 DATA 96x96.xlsx>
>>>>>>
>>>>>
>>>>
>>>
>>
>
NOTICE AND DISCLAIMER
This e-mail (including any attachments) is intended for ...{{dropped:16}}
More information about the Bioconductor
mailing list