Weights for Univariates

src/AverageShiftedHistograms.jl

@@ -2,7 +2,7 @@ module AverageShiftedHistograms
 
 import UnicodePlots, RecipesBase
 using LinearAlgebra, Statistics
-import StatsBase: nobs, fweights
+import StatsBase: nobs, fweights, AbstractWeights
 export ash, ash!, extendrange, xy, xyz, nout, nobs, Kernels
 
 

src/univariate.jl

@@ -23,6 +23,7 @@ function Base.show(io::IO, ::MIME"text/plain", o::Ash)
 end
 
 Base.push!(o::Ash, y::Real) = _histogram!(o::Ash, [y])
+Base.push!(o::Ash, y::Real, weight::Real) = _weightedhistogram!(o::Ash, [y], [weight])
 
 # add data to the histogram
 function _histogram!(o::Ash, y)
@@ -40,6 +41,23 @@ function _histogram!(o::Ash, y)
     return
 end
 
+# add data to the weighted histogram
+function _weightedhistogram!(o::Ash, y, weight::AbstractWeights)
+    b = length(o.rng)
+    a = first(o.rng)
+    b == length(weight) || throw(DimensionMismatch("Length of weights should be same as the length of the range"))
+    δinv = inv(step(o.rng))
+    c = o.counts
+    map(y) do x
+        ki = floor(Int, (x - a) * δinv + 1.5)
+        if 1 <= ki <= b
+            c[ki] += weight[ki]
+        end
+    end
+    o.nobs += length(y)
+    return
+end
+
 # recalculate the ash density
 function _ash!(o::Ash)
     b = length(o.rng)
@@ -62,7 +80,13 @@ end
 # Univariate Ash
     ash(x; kw...)
 
-Fit an average shifted histogram to data `x`.  Keyword options are:
+Fit an average shifted histogram to data `x`.
+
+    ash(x, weight::AbstractWeights; kw...)
+
+Fit a weighted average shifted histogram to data `x`.
+
+Keyword options are:
 
 - `rng`    : values over which the density will be estimated
 - `m`      : Number of adjacent histograms to smooth over
@@ -80,11 +104,12 @@ Fit a bivariate averaged shifted histogram to data vectors `x` and `y`.  Keyword
 - `kernely` : kernel in y direction
 
 # Mutating an Ash object
-Ash objectes can be updated with new data, smoothing parameter(s), or kernel(s).  They cannot, however, change the ranges over which the density is estimated.  It is therefore suggested to err on the side of caution when choosing data endpoints.
+Ash objectes can be updated with new data, new weights(only for univariates), smoothing parameter(s), or kernel(s).  They cannot, however, change the ranges over which the density is estimated.  It is therefore suggested to err on the side of caution when choosing data endpoints.
 
     # univariate
     ash!(obj; kw...)
-    ash!(obj, newx, kw...)
+    ash!(obj, newx; kw...)
+    ash!(obj, newx, weight::AbstractWeights; kw...)
 
     # bivariate
     ash!(obj; kw...)
@@ -97,12 +122,19 @@ function ash(x; nbin=500, rng::AbstractRange = extendrange(x, nbin), m = ceil(In
     _ash!(o)
 end
 
+function ash(x, weight::AbstractWeights; nbin=500, rng::AbstractRange = extendrange(x, nbin), m = ceil(Int, length(rng)/100), kernel = Kernels.biweight)
+    o = Ash(rng, kernel, m)
+    _weightedhistogram!(o, x, weight)
+    _ash!(o)
+end
+
 
 """
     ash!(o::Ash; kw...)
     ash!(o::Ash, newdata; kw...)
+    ash!(o::Ash, newdata, weight::AbstractWeights; kw...)
 
-Update an Ash estimate with new data, smoothing parameter (keyword `m`), or kernel (keyword `kernel`):
+Update an Ash estimate with new data, new weight, smoothing parameter (keyword `m`), or kernel (keyword `kernel`):
 """
 function ash!(o::Ash; m = o.m, kernel = o.kernel)
     o.m = m
@@ -115,6 +147,12 @@ function ash!(o::Ash, y; m = o.m, kernel = o.kernel)
     _histogram!(o, y)
     _ash!(o)
 end
+function ash!(o::Ash, y, weight::AbstractWeights;  m = o.m, kernel = o.kernel)
+    o.m = m
+    o.kernel = kernel
+    _weightedhistogram!(o, y, weight)
+    _ash!(o)
+end
 
 function Base.merge!(o::Ash, o2::Ash)
     o.kernel == o2.kernel || error("Merge failed.  Ash objects use different kernels.")

test/runtests.jl

@@ -107,4 +107,24 @@ end
     xyz(o)
 end
 
+@testset "AshWeighted" begin
+    weight_funcs = (weights, aweights, fweights, pweights)
+
+    for f in weight_funcs
+        x = randn(10_000)
+        o = ash(x, f(ones(21)), rng = -1:0.1:1)
+        o2 = ash(x; rng = -1:0.1:1)
+        @test o == o2
+
+        y = rand(1000)
+        w = f(rand(1:10, 11))
+        o = ash(y; rng = 0:0.1:1)
+        ow = ash(y, w; rng = 0:0.1:1)
+        @test o.counts .* w == ow.counts
+
+        w = f(rand(1:10, 10))
+        @test_throws DimensionMismatch ash(y, w; rng = 0:0.1:1)
+    end 
+end
+
 end #module