[R-sig-Geo] raster stack/brick cor.test p-values

Pascal Oettli kridox at ymail.com
Thu May 29 10:56:33 CEST 2014


Hello Hakim,

The function "cor" can take matrix as input. Not the function "cor.test".

Hope this helps,
Pascal


On Thu, May 29, 2014 at 5:45 PM, Hakim Abdi <hakim.abdi at nateko.lu.se> wrote:
> Dear list:
>
> I am trying to get pixel-wise correlations and significance between two
> raster bricks; one representing monthly NDVI means and another representing
> monthly rainfall means.
>
> I found the following two codes (both from Robert Hijmans) from a previous
> discussions on this list:
> http://r-sig-geo.2731867.n2.nabble.com/questions-on-RasterStack-Brick-td5553580.html
> and
> from StackOverflow:
> http://stackoverflow.com/questions/16698048/r-calculating-pearson-correlation-coefficient-in-each-cell-along-time-line
>
>
> # code 1
>> funcal <- function(xy) {
>>     xy <- na.omit(matrix(xy, ncol=2))
>>     if (ncol(xy) < 2) {
>>         NA
>>     } else {
>>         cor(xy[, 1], xy[, 2])
>>     }
>> }
>> s <- stack(sa, sb)
>> calc(s, funcal)
>> # code 2
>> stackcor <- function(s1, s2, method='spearman') {
>>         mycor <- function(v) {
>>                 x <- v[1:split]
>>                 y <- v[(split+1):(2*split)]
>>                 cor(x, y, method=method)
>>         }
>>         s <- stack(s1, s2)
>>         split <- nlayers(s)/2
>>         calc(s, fun=mycor )
>> }
>
>
> Both codes work as expected with the "cor" function, producing correlation
> grids. However, I tried substituting the cor with cor.test in order to
> extract the p-values:
>
> # using the first code
>> funcal <- function(xy) {
>>     xy <- na.omit(matrix(xy, ncol=2))
>>     if (ncol(xy) < 2) {
>>         NA
>>     } else {
>>         cor.test(xy[, 1], xy[, 2],method="pearson")$p.value
>>     }
>> }
>> s <- stack(sa, sb)
>> calc(s, funcal)
>
>
> I am met with the following error (with trackback in RStudio):
>
> Error in cor.test.default(xy[, 1], xy[, 2], method = "pearson") :
>>   not enough finite observations
>> 8 stop("not enough finite observations")
>> 7 cor.test.default(xy[, 1], xy[, 2], method = "pearson")
>> 6 cor.test(xy[, 1], xy[, 2], method = "pearson")
>> 5 FUN(newX[, i], ...)
>> 4 apply(v, 1, fun)
>> 3 .local(x, fun, ...)
>> 2 calc(meistack, brick.cor.pval)
>> 1 calc(meistack, brick.cor.pval)
>
>
> In the previous r-sig-geo, a user had asked about this error but it was
> left unanswered. I was wondering if anyone can help with this?
>
> Thanks!
>
> Hakim
> ____________________________________________________
>
> Hakim Abdi |  PhD Candidate
> Geobiosphere Science Center
>
> Department of Physical Geography and Ecosystem Science
> Lund University, Sölvegatan 12, SE-223 62 Lund, Sweden
>
> Email: hakim.abdi at nateko.lu.se
>
>         [[alternative HTML version deleted]]
>
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>



-- 
Pascal Oettli
Project Scientist
JAMSTEC
Yokohama, Japan



More information about the R-sig-Geo mailing list