计算相关图的P值

问题描述

我已经将澳大利亚降雨的相关性映射到不同的驱动变量。 stk_b7-是1925-1946年的每月降雨量值,以5公里为单位。这是具有以下尺寸的砖砌堆栈-

stk_b7 
class: RasterStack 
dimensions : 681,841,572721,264  (nrow,ncol,ncell,nlayers)
resolution : 0.05,0.05  (x,y)

然后我有一个csv文件mv_b7,其中包含我的所有驱动程序变量。 SOI,NINO3,NINO4等 列是变量,行是相应时间段的月份。 1925年1月,1925年2月,...

我已使用以下代码计算相关性并将其映射。

myfun <- function(x) cor(as.vector(x),mv_b7$SOI)

#use the function with calc

r <- calc(stk_b7,myfun)

#map results of function applied to raster brick
rasterVis::gplot(r) + 
  geom_tile(aes(fill = value)) + 
   scale_fill_gradient2(low="blue",high="red",mid="white") +
  
  coord_equal() +
  theme_bw() +
  labs(x = "Longitude",y = "Latitude",fill = "Correlation(pearson)",title = "1900-2013 IPO-Precipitation Correlation",subtitle = "COBE sst data provided by NOAA,SILO gridded rainfall") +
  #guides(fill = guide_legend(override.aes = list(size = 2),ncol = 2)) +
  theme(plot.title = element_text(size = 12,hjust = 0.5,face = "bold"),plot.subtitle = element_text(size = 8,hjust = 0.5),legend.title = element_text(size = 7),legend.text = element_text(size = 7),axis.title = element_text(size = 10),axis.text = element_text(size = 7))

这将根据需要生成关联图。但是,现在我需要计算相应的P值。我不确定这是我要映射的东西还是整个关联图只得到一个P值。我尝试使用具有相同格式的rcorr和cor.test,但收到以下错误

    myfun <- function(x) cor.test(as.vector(x),mv_b7$SOI)
    r <- calc(stk_b7,myfun)
    Error in .calcTest(x[1:5],fun,na.rm,forcefun,forceapply) : 
    cannot use this function

不确定这是否与NA有关,但是对于为什么正常相关将起作用而其他方面却不知道却感到困惑。对此的任何帮助将不胜感激

解决方法

这是因为cor.test的输出具有多个值,并且需要明确使用cor.test才能获得该函数调用的p值。

a <- rnorm(4) #Generating fake data
b <- rnorm(4)

cor.test(a,b)

Pearson's product-moment correlation

data:  a and b
t = -0.47545,df = 2,p-value = 0.6813
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
 -0.9797036  0.9260328
sample estimates:
       cor 
-0.3186697 

要实际访问p值,您可以在cor.test的功能文档中看到,您将需要访问输出中的第三个值。

cor.test(a,b)[3]
$p.value
[1] 0.6813303

因此,我将更改您的函数以确保您在cor.test()调用中访问第三个值,然后它应该可以工作。