[R] help with function

[Iain Gallagher]

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.

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
>