CODEWITHSUNDEEP
statisticsjuliamoving-average
clearseplex
clearseplex

Reputation: 719

Moving average in Julia

I want to calculate the simple moving average of a array in Julia. I have a plain simple array, but all packages I found require a TimeArray to calculate the moving average. Is there a package that doesn't require me to artificially create a TimeArray?

Upvotes: 5

Views: 8331

Answers (6)

Ammar Bandukwala
Ammar Bandukwala

Reputation: 1582

Here's an implementation using Base.foldl:

function sma(in::Vector{<:Number}, windowLen::Integer)::Vector{Number}
  window = []
  return foldl(
    (
      (acc, v) -> begin
        push!(window, v)
        if size(window, 1) > windowLen
          popfirst!(window)
        end
        return push!(acc, mean(window))
      end
    ), in, init=[]
  )
end

Example result:

sma([1, 2, 3, 4, 5], 2) = Number[1.0, 1.5, 2.5, 3.5, 4.5]

Upvotes: 1

Max
Max

Reputation: 835

Here is a two-liner that becomes a one-liner in the (likely) event that you already have a mean function somewhere.

using Statistics: mean
movingaverage(g, n) = [i < n ? mean(g[begin:i]) : mean(g[i-n+1:i]) for i in 1:length(g)]

Like some of the longer answers posted, this has the advantage of preserving the length of the input vector g, which is useful for plotting.

Upvotes: 1

Satvik Beri
Satvik Beri

Reputation: 271

I tried implementing a few versions:

function rolling_sum(arr, n)
    so_far = sum(arr[1:n])
    out = zero(arr[n:end])
    out[1] = so_far
    for (i, (start, stop)) in enumerate(zip(arr, arr[n+1:end]))
        so_far += stop - start
        out[i+1] = so_far
    end
    return out
end

rolling_mean(arr, n) = rolling_sum(arr, n) ./ n

function rolling_mean2(arr, n)
    return imfilter(arr, OffsetArray(fill(1/n, n), -n), Inner())
end

function rolling_mean3(arr, n)
    so_far = sum(arr[1:n])
    out = zero(arr[n:end])
    out[1] = so_far
    for (i, (start, stop)) in enumerate(zip(arr, arr[n+1:end]))
        so_far += stop - start
        out[i+1] = so_far / n
    end
    return out
end

function rolling_mean4(arr, n)
    rs = cumsum(arr)[n:end] .- cumsum([0.0; arr])[1:end-n]
    return rs ./ n
end

julia> v = rand(50_000);
julia> @btime(rolling_mean($v, 200));
  125.362 μs (9 allocations: 1.52 MiB)

julia> @btime(rolling_mean2($v, 200));
  1.089 ms (127 allocations: 3.06 MiB)

julia> @btime(rolling_mean3($v, 200));
  93.137 μs (7 allocations: 1.14 MiB)

julia> @btime(rolling_mean4($v, 200));
  161.613 μs (12 allocations: 2.28 MiB)

In this case, it looks like the specialized version that adds and removes an element at each step can be significantly faster.

Upvotes: 2

tholy
tholy

Reputation: 12179

julia> using ImageFiltering, OffsetArrays

julia> v = zeros(20); v[10] = 1   # test vector to perform moving average
1

julia> kernel = OffsetArray(fill(1/8, 8), -5:2)  # moving average of 5 previous, current, and 2 ahead
8-element OffsetArray(::Array{Float64,1}, -5:2) with eltype Float64 with indices -5:2:
 0.125
 0.125
 0.125
 0.125
 0.125
 0.125
 0.125
 0.125

julia> [v imfilter(v, kernel)]
20×2 Array{Float64,2}:
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.125
 0.0  0.125
 1.0  0.125
 0.0  0.125
 0.0  0.125
 0.0  0.125
 0.0  0.125
 0.0  0.125
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.0  
 0.0  0.0

Upvotes: 4

Steven Siew
Steven Siew

Reputation: 873

you can write your own moving average like I did

function movingaverage(X::Vector,numofele::Int)
    BackDelta = div(numofele,2) 
    ForwardDelta = isodd(numofele) ? div(numofele,2) : div(numofele,2) - 1
    len = length(X)
    Y = similar(X)
    for n = 1:len
        lo = max(1,n - BackDelta)
        hi = min(len,n + ForwardDelta)
        Y[n] = mean(X[lo:hi])
    end
    return Y
end

Upvotes: 2

Przemyslaw Szufel
Przemyslaw Szufel

Reputation: 42194

What about:

moving_average(vs,n) = [sum(@view vs[i:(i+n-1)])/n for i in 1:(length(vs)-(n-1))]

This could be further optimized by making a standard for loop, pre-allocating the result array and at each iteration subtracting and adding just one element of the input array. However, for most applications the above simple code is sufficient.

Upvotes: 10

Related Questions