数组中每12个矩阵的逐元素平均值,重复一个序列12次,无需使用for循环。

3

我有一个维度为[360, 180, 396]的数组。这些是经度、纬度和月份-年份,用于33年的月度数据。元素是该纬度/经度的百分比。

从这个数组中,我想制作一个摘要数组,以便在后续分析中使用,而不是默认使用for循环。我想得到每个月所有33年数据的平均值,然后是所有年份的年平均值。 这是我制作的包含数据的摘要数组的空白。

mca <- array(data = NA, 
             dim = c(360,180,13), 
             dimnames = list(lon, 
                             lat, 
                             c(month.abb, "Ann")))

以下是这个例子的较小测试输入和输出数组

#input

set.seed(42)
smallin <- array(data = rnorm(n = 600, mean = 60, sd = 20),
               dim = c(5, 5, 24))


#output to fill
smallout <- array(data = NA, 
             dim = c(5,5,13), 
             dimnames = list(c("1", "2", "3", "4", "5"), 
                             c("-89.5", "-88.5", "-87.5", "-86.5", "-85.5"), 
                             c(month.abb, "Ann")))

根据这个问题的第二个答案,我尝试了

jan <- apply(ca, c(seq(from = 1, to = 385, by = 12)), mean)

#also 

ind_jan <- c(seq(from = 1, to = 385, by = 12))
jan <- apply(ca, ind_jan, mean)

我认为这相当于

jan <- apply(smallin, c(seq(from = 1, to = 13, by = 12)), mean)

考虑到边距,我需要放置第三个维度以进行平均,但是出现了错误:

apply(ca, c(seq(from = 1, to = 385, by = 12)), mean) 中的错误: 'MARGIN' does not match dim(X)

我回到上面的查询并意识到 margin = 1:2 必须选择每个矩阵的所有部分(尺寸为 1 和 2)。因此,使用它,我可以得到所有矩阵的平均值,这应该是输出数组 [,,13] 的年平均百分比。

smallout[,,13] <- apply(smallin, 1:2, mean)

但我仍然不知道如何让它从第1个矩阵开始,每12个矩阵取平均值,然后从第2个矩阵开始,再从第3个矩阵开始...

我已经阅读了apply文档,但在这种情况下/难以理解。所有建议的问题似乎都是用Python(或其他语言)。

我也不确定是否可以一次完成所有操作,还是通过索引将矩阵逐个传递到输出数组中。

我能想到的最接近的方法是类似于

ind_jan <- c(seq(from = 1, to = 13, by = 12))
smallout[,,1] <- apply(smallin[,,c(ind_jan)], 1:2, mean)

对于数组中的每个输出矩阵重复。是否有更少手动/更高效/更好的方法?


你能否提供一个的可重现示例,比如3x4x2而不是360x180x13,并使用示例数据填充它,而不是NA,然后展示该示例输入的期望输出?我认为这将使您和我们更容易理解正在发生的事情,并为我们提供可测试的东西。 - Gregor Thomas
另外,如果不知道“lon”和“lat”,我们也无法运行您的代码。 - jay.sf
我已经添加了一个迷你示例,并提供了可能/部分解决方案(仍不是非常高效),这是我在向别人解释时想到的。 - The_Tams
3个回答

3

考虑这个简化的数组A(见下面的数据)。

str(A)
# int [1:2, 1:3, 1:6] 1 1 1 1 1 1 2 2 2 2 ...

我们可以使用 sapply 来“循环”处理每一年,同时选项 simplify='array' 可以返回一个年度平均数数组。
yrs <- seq_len(dim(A)[3]/nm)
sapply(yrs, \(i) apply(A[, , 1:nm + i - 1], 1:2, mean), simplify='array')
# , , 1
# 
#      [,1] [,2] [,3]
# [1,]    2    2    2
# [2,]    2    2    2
# 
# , , 2
# 
#      [,1] [,2] [,3]
# [1,]    2    2    2
# [2,]    2    2    2

因此,相应地,跨年度的月平均数:

mnt <- seq_len(nm)
sapply(mnt, \(i) apply(A[, , i], 1:2, mean), simplify='array')
# , , 1
# 
#      [,1] [,2] [,3]
# [1,]    1    1    1
# [2,]    1    1    1
# 
# , , 2
# 
#      [,1] [,2] [,3]
# [1,]    2    2    2
# [2,]    2    2    2
# 
# , , 3
# 
#      [,1] [,2] [,3]
# [1,]    3    3    3
# [2,]    3    3    3

数据:

nm <- 3  ## no. "months"  ## actually 12 months in real years
ny <- 2  ## no. "years"  ## in your case 33
A <- array(rep(1:nm, each=nm*ny), c(2, 3, nm*ny))  ## think this is your `ca`

当我尝试将此代码应用于我的示例代码(smallin)时,yrs <- seq_len(dim(smallin)[3]/12) sapply(yrs, \(i) apply(smallin[, , i], 1:2, mean), simplify='array')我的输出是两个矩阵,这表明这是每年的平均值而不是每月的平均值。而且,mnt <- seq_len(12) sapply(mnt, \(i) apply(smallin[, , 1:12 + i - 1], 1:2, mean),simplify='array')结果是12个矩阵。所以年度和月份的代码块是否颠倒了?如果我对数组进行切片,我也不确定您所说的失去一个维度的含义是什么... - The_Tams
1
@The_Tams 是的,i 指的是函数 \(i)(在 R>4.1 中是 function(i) 的简写)。 - jay.sf
@The_Tams 可能我们理解直径的术语,例如“年平均数”,我指的是每个月的平均数,所以我们有三个(就像简化后的月份数)。 - jay.sf
好的,你的意思是年平均值是每个月33年数据的平均值,而月平均值是一年中所有月份的平均值。 我想我需要更多地研究sapply函数内部的公式,因为虽然我可以看到它似乎有效,但我不理解它。 - The_Tams
让我们在聊天中继续这个讨论 - The_Tams
显示剩余4条评论

2

您可以通过将包含月份和年份的最后一个维度拆分为单独的月份和年份维度,为数组添加另一个维度。

i <- dim(smallin)
dim(smallin) <- c(i[1:2], 12L, i[3]/12L)

使用这个方法,你可以得到每个月在所有年份中的平均值:
apply(smallin, 1:3, mean)
#, , 1
#
#         [,1]     [,2]     [,3]     [,4]     [,5]
#[1,] 73.66338 58.35988 72.33907 62.19628 52.08766
#[2,] 61.95544 79.93891 75.27725 49.30859 44.07820
#[3,] 64.02119 68.98285 35.76780 35.06961 58.79089
#[4,] 73.67935 67.72028 50.90479 23.22819 72.14434
#[5,] 62.57796 59.03798 64.53486 83.65987 97.04576
#
#...
#
#, , 12
#
#         [,1]     [,2]     [,3]     [,4]     [,5]
#[1,] 83.55254 68.77645 48.88358 52.99573 56.82992
#[2,] 83.47723 39.02472 95.08051 65.97988 54.00097
#[3,] 47.59936 36.93396 38.35189 57.86126 83.99976
#[4,] 73.00906 53.71818 36.93229 80.85843 39.27094
#[5,] 81.67441 64.50031 62.71359 56.27758 54.01388

单年度的年均值:

apply(smallin, c(1,2,4), mean)
#, , 1
#
#         [,1]     [,2]     [,3]     [,4]     [,5]
#[1,] 60.77253 60.15417 54.71206 67.31820 62.05012
#[2,] 56.60298 59.14604 73.17469 57.66912 53.36540
#[3,] 56.52924 56.31096 58.73874 67.47850 59.06819
#[4,] 67.75999 56.45636 49.43743 55.14660 65.46497
#[5,] 60.28056 62.17656 55.08681 54.15788 60.05240
#
#, , 2
#
#         [,1]     [,2]     [,3]     [,4]     [,5]
#[1,] 60.55035 65.21223 59.92112 59.75500 69.77088
#[2,] 60.89782 54.59722 55.17699 59.06815 60.03906
#[3,] 58.85733 54.02893 47.31326 63.10434 59.56569
#[4,] 60.96362 61.82648 55.45109 54.50272 45.21176
#[5,] 59.94452 54.31497 60.64839 64.65777 80.86525

所有年份的年度平均值:

apply(smallin, 1:2, mean)
#         [,1]     [,2]     [,3]     [,4]     [,5]
#[1,] 60.66144 62.68320 57.31659 63.53660 65.91050
#[2,] 58.75040 56.87163 64.17584 58.36864 56.70223
#[3,] 57.69329 55.16994 53.02600 65.29142 59.31694
#[4,] 64.36180 59.14142 52.44426 54.82466 55.33836
#[5,] 60.11254 58.24577 57.86760 59.40782 70.45883

谢谢你,GKi。我之前没意识到改变现有数组的形状是如此容易。我喜欢这个答案让原始数组更直观,方便我处理数据,以及想要做的其他任何事情。 - The_Tams

1

我相信有更好的方法,如果有人知道,请告诉我。但是,通过我研究如何进行索引以选择每个月数据并使用apply取平均值,下面的方法确实有效。

mca[,,1] <- apply(ca[,,c(seq(from = 1, to = 396, by = 12))], 1:2, mean)
mca[,,2] <- apply(ca[,,c(seq(from = 2, to = 396, by = 12))], 1:2, mean)
mca[,,3] <- apply(ca[,,c(seq(from = 3, to = 396, by = 12))], 1:2, mean)
mca[,,4] <- apply(ca[,,c(seq(from = 4, to = 396, by = 12))], 1:2, mean)
mca[,,5] <- apply(ca[,,c(seq(from = 5, to = 396, by = 12))], 1:2, mean)
mca[,,6] <- apply(ca[,,c(seq(from = 6, to = 396, by = 12))], 1:2, mean)
mca[,,7] <- apply(ca[,,c(seq(from = 7, to = 396, by = 12))], 1:2, mean)
mca[,,8] <- apply(ca[,,c(seq(from = 8, to = 396, by = 12))], 1:2, mean)
mca[,,9] <- apply(ca[,,c(seq(from = 9, to = 396, by = 12))], 1:2, mean)
mca[,,10] <- apply(ca[,,c(seq(from = 10, to = 396, by = 12))], 1:2, mean)
mca[,,11] <- apply(ca[,,c(seq(from = 11, to = 396, by = 12))], 1:2, mean)
mca[,,12] <- apply(ca[,,c(seq(from = 12, to = 396, by = 12))], 1:2, mean)
mca[,,13] <- apply(ca, 1:2, mean)

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接