On Fri, 17 Dec 2010, Iain Gallagher wrote:
> Hello List
>
> I'm moving this over from the bioC list as, although the problem
I'm working on is biological, the current bottle neck is my poor
understanding of R.
>
> I wonder if someone would help me with the following function.
>
Here is how I'd take it apart.
Either
1) put browser() as the first line of the function,then feed lines to the
browser one -by-one to see where it hangs,
2) trace(cumulMetric) , then try to run it (much like 1, but it will
handle feeding the lines of the function more easily)
3) options( error = recover ), then run it till it hits the error, then
browser thru the frames to see what is where
See
?browser
?trace
?recover
as background.
HTH,
Chuck
> cumulMetric <- function(deMirPresGenes, deMirs){
> ???
> #need to match position of each miR in deMirPresGenes with its FC to form a
vector of FC in correct order
> ? ? fc <- deMirs
> ? ? fcVector <- as.numeric(with (fc, FC[match(deMirPresGenes[,4],
Probe)] ) )
>
> ? ? #multiply fc by context score for each interaction
> ? ? metric <- fcVector * as.numeric(deMirPresGenes[,11])
> ? ? geneMetric <- cbind(deMirPresGenes[,2], as.numeric(metric))
>
> ? ? #make cumul weighted score
> ? ? listMetric <- unstack(geneMetric,
as.numeric(geneMetric[,2])~geneMetric[,1])
> ? ? listMetric <- as.data.frame(sapply(listMetric,sum)) #returns a
dataframe
> ? ? colnames(listMetric) <- c('cumulMetric')
>
> ? ? #return whole list
> ? ? return(listMetric)
> }
>
> deMirPresGenes looks like this:
>
> Gene.ID? ? Gene.Symbol? ? Species.ID? ? miRNA? ? Site.type? ? UTR_start? ?
UTR_end? ? X3pairing_contr? ? local_AU_contr? ? position_contr? ? context_score?
? context_percentile
> 22848? ? AAK1? ? 9606? ? hsa-miR-183? ? 2? ? 1546? ? 1552? ? -0.026? ?
-0.047? ? 0.099? ? -0.135? ? 47
> 19? ? ABCA1? ? 9606? ? hsa-miR-183? ? 2? ? 1366? ? 1372? ? -0.011? ?
-0.048? ? 0.087? ? -0.133? ? 46
> 20? ? ABCA2? ? 9606? ? hsa-miR-495? ? 2? ? 666? ? 672? ? -0.042? ? -0.092?
? -0.035? ? -0.33? ? 93
> 23456? ? ABCB10? ? 9606? ? hsa-miR-183? ? 3? ? 1475? ? 1481? ? 0.003? ?
-0.109? ? -0.05? ? -0.466? ? 98
> 6059? ? ABCE1? ? 9606? ? hsa-miR-495? ? 2? ? 1474? ? 1480? ? 0.005? ?
-0.046? ? 0.006? ? -0.196? ? 58
> 55324? ? ABCF3? ? 9606? ? hsa-miR-1275? ? 3? ? 90? ? 96? ? 0.007? ? 0.042?
? -0.055? ? -0.316? ? 94
>
> although it is much longer in 'real life'.
>
> The aim of the function is to extract a dataframe of gene symbols along
with a weighted score from the above data. The weighted score is the FC column
of deMirs * the context_score column of deMirPresGenes. This is easy peasy!
>
> Where I'm falling down is that if I run this function it complains that
'geneMetric' can't be found. Hmm - I've run it all line by line
(i.e. not as a function) and it works but wrapped up like this it fails!
>
> e.g.
>
>> testF2 <- cumulMetric(testF1, deMirs$up)
> Error in eval(expr, envir, enclos) : object 'geneMetric' not found
>
> deMirs$up looks like this:
>
> Probe? ? FC
> hsa-miR-183? ? 2.63
> hsa-miR-1275? ? 2.74
> hsa-miR-495? ? 3.13
> hsa-miR-886-3p? ? 3.73
> hsa-miR-886-5p? ? 3.97
> hsa-miR-144*? ? 6.62
> hsa-miR-451? ? 7.94
>
> In an effort to debug this I have examined each object using
'print' statements (as suggested by cstrato on the bioC list).
>
> All the objects in the function up until listMetric look ok in terms of
structure and length. But the error persists and listMetric is not made?!?! Odd.
>
> I have added some comments to the output below.
>
>> tf2<-cumulMetric(tf1, deMirs$up)#deMirs$up is a dataframe (see above
- it is the same as deMirs)
>
>
> [1] 2.63 2.63 3.13 2.63 3.13 2.74 # print fcVector - looks ok
> [1] -0.35505 -0.34979 -1.03290 -1.22558 -0.61348 -0.86584 # print metric -
looks ok
> [1] 1045 # length of metric - is correct
> ? ???sym? ? ? metric? ? # print geneMetric - looks ok
> [1,] "AAK1"???"-0.35505"
> [2,] "ABCA1"? "-0.34979"
> [3,] "ABCA2"? "-1.0329"
> [4,] "ABCB10" "-1.22558"
> [5,] "ABCE1"? "-0.61348"
> [6,] "ABCF3"? "-0.86584"
> [1] 1045 # nrow of geneMetric - is correct
> Error in eval(expr, envir, enclos) : object 'geneMetric' not found
>>
>
> Could someone possibly point out where I falling down.
>
> Thanks
>
> iain
>
>> sessionInfo()
> R version 2.12.0 (2010-10-15)
> Platform: x86_64-pc-linux-gnu (64-bit)
>
> locale:
> [1] LC_CTYPE=en_GB.utf8? ? ???LC_NUMERIC=C? ? ? ? ? ???
> [3] LC_TIME=en_GB.utf8? ? ? ? LC_COLLATE=en_GB.utf8???
> [5] LC_MONETARY=C? ? ? ? ? ???LC_MESSAGES=en_GB.utf8???
> [7] LC_PAPER=en_GB.utf8? ? ???LC_NAME=C? ? ? ? ? ? ???
> [9] LC_ADDRESS=C? ? ? ? ? ? ? LC_TELEPHONE=C? ? ? ? ???
> [11] LC_MEASUREMENT=en_GB.utf8 LC_IDENTIFICATION=C? ???
>
> attached base packages:
> [1] stats? ???graphics? grDevices utils? ???datasets? methods???base? ???
>
> loaded via a namespace (and not attached):
> [1] tools_2.12.0
>>
>
>
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>
Charles C. Berry Dept of Family/Preventive Medicine
cberry at tajo.ucsd.edu UC San Diego
http://famprevmed.ucsd.edu/faculty/cberry/ La Jolla, San Diego 92093-0901