[kth-package] Re: how to do with sevral replicates?
- From: Valtteri Wirta <valtteri@xxxxxxxxxxxxxx>
- To: kth-package@xxxxxxxxxxxxx
- Date: Fri, 18 Nov 2005 14:06:29 +0100
Hej Monica,
Below is the code for calcMeanSD2 function that probably does what you need.
I wrote and tested it on the gpr file you send
and it seems to work nicely. You must however
test it and manually check a couple of the
calculations... Also note that function does not
weight probes differently if they are on the same array or on different arrays.
The arguments that need to be defined in the function call are:
ma.obj = the ma-object containing the M-values
slideVec = a numeric vector that contains the
slide numbers that you wish to include in the calculations
ids = a character vector containing the names of
all probes on the array. Must be in the same
order as the M-values in the ma.obj. Use for example gpr$Name[,1] here
calcMeanSD2 <- function(ma.obj, slideVec, ids)
{
# Create a list of unique ids in the input list
uniq.ids <- unique(ids)
# Concatenate all M-values and ids into long
list. Include only data from slides that are defined in the slideVec
Mvec <- vector(mode="numeric")
Namevec <- vector(mode="character")
for(counter in slideVec)
{
Namevec <- c(Namevec, ids)
Mvec <- c(Mvec, ma.obj$M[,counter])
}
# for each uniq.ids, identify its position
(index) in the long list that was created above
M.locs <- list()
for(counter2 in 1:length(uniq.ids))
{
M.locs[[counter2]] <- which(Namevec == uniq.ids[counter2])
}
# calculate the mean, sd and n for each uniq.ids
och and enter these into a matrix containing nrow=lenght(uniq.ids), ncol=3
out <- matrix(ncol=3, nrow=length(uniq.ids), data=NA)
colnames(out) <- c("mean M", "SD M", "count")
rownames(out) <- uniq.ids
for(counter3 in 1:length(uniq.ids))
{
out[counter3, "mean M"] <- mean(Mvec[M.locs[[counter3]]], na.rm=TRUE)
out[counter3, "SD M"] <- sd(Mvec[M.locs[[counter3]]], na.rm=TRUE)
out[counter3, "count"] <-
length(which(is.na(Mvec[M.locs[[counter3]]])==FALSE))
}
# summarise the results and return the out object
windows()
subplots(4)
hist(out[,"mean M"], nclass=50, col=4, main="M-values")
hist(out[,"SD M"], nclass=50, col=4, main="SD")
hist(out[,"count"], col=4, main="Counts", nclass=100)
count.nums <- vector(mode="numeric")
for(counter4 in 1:10)
{
count.nums[counter4] <- length(which(out[,"count"] == counter4))
}
barplot(count.nums, col=4, names=1:10, main="Counts")
return(out)
}
Example of usage
results.out <- calcMeanSD2(ma.obj=ma, slideVec=1:3, ids=gpr$Name[,1])
good luck,
VW
At 15:22 2005-11-17, you wrote:
> Hej Valtteri!
Vi har kollat på skriptet "reorderMalaria" som du och Uffe gjorde. Men det
är ändå ganska komplicerat att få det helt rätt. Det är ju avgörande att
utmatriserna blir helt rätt. Vi har tre replikat överallt, dvs inga
multipler men vi har många "empty" och "blank" som inte ska ingå i
analysen. Vi skulle alltså behöva sortera om ordningen så att R räknar på
triplikat. Jag har suttit ganska länge och försökt skriva om calcMeanSD
funktionen men som tur var så pratade jag med Uffe idag som sa just att
man måste sortera om dem. Jag pluggar just nu på KTH och har läst lite
java men det är ändå ganska krångligt. Vi vill kunna beräkna medelvärdet
av M-värdena och SD för tre replikat för alla slides inom varje stam. Dvs
den ska retunera om man tex har fyra slides för en stam ett hopslaget
medelvärde av M-värdena per genID samt ett hopslaget SD värde för varje
gen. Vi har oftast tre slides per stam några fyra. Jag bifogar en gpr-fil
av stammen SME33 som ska ingå i analysen så ser du hur den är upplagd. Jag
uppskattar verkligen om du har tid och titta på detta.
Uffe skickar också en puss i ljumsken!
Mvh
Monica A, Christel Blomberg, Bakt 3, SMI
Hi Monica,
>
> calcMeanSD can't probably handle that.
>
> There are many ways of doing this, however they may require some
> knowledge of the way R does the indexing...
>
> If your replicate probes are in a regular order (same number of rows
> between each) then this is an easy calculation.
> If they are not, then you may have to somehow re-order the probes and
> then do the calculations... if you send me an example file off-list I
> can have a look at it.
>
> regards,
> Valtteri
>
>
>
> At 22:59 2005-11-16, you wrote:
>>Hi, I wonder if you could use the R-function calcMeanSD for more than
>>duplicates. I've tried to change the code but still the function only
>>count for duplicates. Do it exist some other function that calculate the
>>same but for sevral replicates?
>>
>>H/Monica
>>
>>**********************************************
******************************
>>Mailing list for KTH-package. For more information see:
>>www.biotech.kth.se/molbio/microarray
>>To change the settings or unsubscribe:
>>//www.freelists.org/cgi-bin/list?list_id=kth-package
>>To send mail to the list use: kth-package@xxxxxxxxxxxxx
>>To get information on how the mailing list works, send an empty
>>e-mail with the word 'faq' to kth-package-request@xxxxxxxxxxxxx
>>***************************************************************************
>
>
>
****************************************************************************
> Mailing list for KTH-package. For more information see:
> www.biotech.kth.se/molbio/microarray
> To change the settings or unsubscribe:
> //www.freelists.org/cgi-bin/list?list_id=kth-package
> To send mail to the list use: kth-package@xxxxxxxxxxxxx
> To get information on how the mailing list works, send an empty e-mail
> with the word 'faq' to kth-package-request@xxxxxxxxxxxxx
> ***************************************************************************
>
>
****************************************************************************
Mailing list for KTH-package. For more information see: www.biotech.kth.se/molbio/microarray
To change the settings or unsubscribe: //www.freelists.org/cgi-bin/list?list_id=kth-package
To send mail to the list use: kth-package@xxxxxxxxxxxxx
To get information on how the mailing list works, send an empty e-mail with the word 'faq' to kth-package-request@xxxxxxxxxxxxx
***************************************************************************
Other related posts:
