Julia 中矩阵列的平均值

Average over the columns of the matrix in Julia

我有一个大矩阵,其中包含以下形式的浮点条目

[ a b c d 
  e f g h
  i j k l
  m n o p ]

一些值是离群值,所以我想用相应列中最近的 k 个条目对每个条目取平均值并保留形状。换句话说,对于 k = 3:

有这样的东西
[        a                  b                   c                  d 
      (e + a)/2          (f + b)/2          (g + c)/2          (h + d)/2
    (e + a + i)/3      (f + b + j)/3      (g + c + k)/3      (h + d + l)/3
    (e + i + m)/3      (f + j + n)/3      (g + k + o)/3      (h + l + p)/3   ] 


 etc.

您可以使用 RollingFunctionsmapslices 来做到这一点:

julia> a = reshape(1:16, 4, 4)
4×4 reshape(::UnitRange{Int64}, 4, 4) with eltype Int64:
 1  5   9  13
 2  6  10  14
 3  7  11  15
 4  8  12  16

julia> using RollingFunctions

julia> mapslices(x -> runmean(x, 3), a, dims = 1)
4×4 Matrix{Float64}:
 1.0  5.0   9.0  13.0
 1.5  5.5   9.5  13.5
 2.0  6.0  10.0  14.0
 3.0  7.0  11.0  15.0

我不知道 RollingFunctions,但常规循环要快 4 倍。不知道是不是mapslices导致的某种类型不稳定?

function runmean(a,W)
    A = similar(a)
    for j in axes(A,2), i in axes(A,1)
        l = max(1, i-W+1)
        A[i,j] = mean(a[k,j] for k=l:i)
    end
    A 
end

测试产量:

using RollingFunctions 
@btime mapslices(x -> runmean(x, 3), A, dims = 1) setup=(A = rand(0.0:9,1000,1000))
@btime runmean(A,3) setup=(A = rand(0.0:9,1000,1000))

  15.326 ms (10498 allocations: 23.45 MiB)
   4.410 ms (2 allocations: 7.63 MiB)