Added mase error measure for continuous variables (Mean Absolute Scaled Error)

Also:
- improved perfomances of the Scaler algorithm using inplace assignments
- perform random reset of a NN model that doesn't converge  only if the model has not beeing fitted before

Author: sylvaticus

src/Nn/Nn.jl

@@ -1178,9 +1178,9 @@ function fit!(m::NeuralNetworkEstimator,X,Y)
         (onfail == "stop") && break # break as one training done, whatever the result
         if a == fail_attempts
             error("Fitting the model on the data failed. Loss is not decreasing even after $(fail_attempts) attempts.")
-        else
+        elseif ! m.fitted # random reset only if the model has not beeing fitted before
             layers = m.par.nnstruct.layers
-            random_init!.(layers;rng=rng)
+            random_init!.(layers;rng=rng) 
         end    
     end
 

src/Utils/Measures.jl

@@ -749,6 +749,50 @@ function relative_mean_error(y,ŷ;normdim=false,normrec=false,p=1)
 end
 
 
+"""
+  mase(y, ŷ;seq=false)
+
+Mean absolute scaled error between y and ŷ.
+
+A relative measure of error that is scale invariant.
+
+If `seq` is `true`, the scaling is done using the in-sample one-step naive forecast method, and the function returns the Mean Absolute Error (MAE) divided by the mean absolute deviation of the difference between consecutive observations of y.
+
+When `seq` is `false` (default), the scaling is done using the in-sample mean and the function returns the Mean Absolute Error (MAE) divided by the mean absolute deviation of y. This is also equal to the relative mean error of the zero-mean, unit variance scaled data.
+
+See [Wikipedia article](https://en.wikipedia.org/wiki/Mean_absolute_scaled_error).
+
+# Notes
+- this function works only for vector data (i.e. single output regression)
+- use `seq=true` for sequential data (e.g. time series)
+
+# Example
+```julia
+julia> using BetaML
+julia> y = [-2, 5, 6, 0, 3, 5, -10, 15, -5];
+julia> ŷ = [0, 4, 5, 1, 2, 4, -8, 18, -6];
+julia> mean_absolute_scaled_error = BetaML.mase(y,ŷ)
+0.2647058823529412
+julia> ms = BetaML.Scaler()
+A Scaler BetaMLModel (unfitted)
+julia> y_s = BetaML.fit!(ms, y);
+julia> ŷ_s = BetaML.predict(ms, ŷ);
+julia> mean_absolute_scaled_error_s = BetaML.relative_mean_error((y_s,ŷ_s)
+0.2647058823529412
+```
+"""
+function mase(y::AbstractVector,ŷ::AbstractVector;seq=false)
+    if seq
+        n = length(y)
+        d = sum(abs.(y[2:end] .- y[1:end-1]))/(n-1)
+    return sum(abs.(y .- ŷ))/ (n*d)
+    else
+        ȳ = mean(y)
+        d = sum(abs.(y .- ȳ))
+        return sum(abs.(y .- ŷ)) / d
+    end
+end
+
 # ------------------------------------------------------------------------------
 # FeatureRanker
 

src/Utils/Processing.jl

@@ -646,7 +646,7 @@ function _fit(m::MinMaxScaler,skip,X,cache)
           imin,imax =   (m.inputRange[1](skipmissing(c)),  m.inputRange[2](skipmissing(c)) )
           if cache
             omin, omax = m.outputRange[1], m.outputRange[2]
-            X_scaled[:,ic] = (c .- imin) .* ((omax-omin)/(imax-imin)) .+ omin
+            X_scaled[:,ic] .= (c .- imin) .* ((omax-omin)/(imax-imin)) .+ omin
           end
           push!(actualRanges,(imin,imax))
         else
@@ -675,7 +675,7 @@ function _fit(m::StandardScaler,skip,X::AbstractArray,cache)
             sfμ[ic] = - μ
             sfσ[ic] = 1 ./ sqrt.(σ²)
             if cache
-                X_scaled[:,ic] = (c .+ sfμ[ic]) .* sfσ[ic] 
+                X_scaled[:,ic] .= (c .+ sfμ[ic]) .* sfσ[ic] 
             end
         end
     end
@@ -715,15 +715,15 @@ function _predict(m::StandardScaler,pars::StandardScaler_lp,skip,X;inverse=false
         xnew .= X
         for (ic,c) in enumerate(eachcol(X))
             if !(ic in skip)
-                xnew[:,ic] = (c .+ pars.sfμ[ic]) .* pars.sfσ[ic] 
+                xnew[:,ic] .= (c .+ pars.sfμ[ic]) .* pars.sfσ[ic] 
             end
         end
         return xnew
     else
         xorig = deepcopy(X)
         for (ic,c) in enumerate(eachcol(X))
             if !(ic in skip)
-                xorig[:,ic] = (c ./ pars.sfσ[ic] .- pars.sfμ[ic])  
+                xorig[:,ic] .= (c ./ pars.sfσ[ic] .- pars.sfμ[ic])  
             end
         end
         return xorig

src/Utils/Utils.jl

@@ -52,7 +52,7 @@ export @codelocation, generate_parallel_rngs,
        dtanh, sigmoid, dsigmoid, softmax, dsoftmax, dmaximum, dmean, pool1d, softplus, dsoftplus, mish, dmish, # exp/trig based functions
        bic, aic,
        autojacobian, match_known_derivatives,
-       squared_cost, dsquared_cost, mse, crossentropy, dcrossentropy, class_counts, kl_divergence, class_counts_with_labels, mean_dicts, mode, gini, entropy, variance, sobol_index,
+       squared_cost, dsquared_cost, mse, mase, crossentropy, dcrossentropy, class_counts, kl_divergence, class_counts_with_labels, mean_dicts, mode, gini, entropy, variance, sobol_index,
        error, accuracy, relative_mean_error,
        ConfusionMatrix, ConfusionMatrix_hp, silhouette,
        cross_validation,

test/Utils_tests.jl

@@ -290,6 +290,62 @@ avgϵRel = (sum(abs.(ŷ-y).^p)^(1/p) / (n*d)) / (sum( abs.(y) .^p)^(1/p) / (n*d
 @test relative_mean_error(y,ŷ,normdim=false,normrec=false,p=p) == avgϵRel
 
 
+@testset "mase (Mean Absolute Scaled Error) test" begin
+    obs = [-2, 5, 6, 0, 3, 5, -10, 15, -5]
+    est = [0, 4, 5, 1, 2, 4, -8, 18, -6]
+    obs2 = 100 .+ obs / 10
+    est2 = 100 .+ est / 10
+    obsb = vcat(obs,obs,obs,obs)
+    estb = vcat(est,est,est,est)
+    obs2b = vcat(obs2,obs2,obs2,obs2)
+    est2b = vcat(est2,est2,est2,est2)
+
+    rme = BetaML.relative_mean_error(obs,est)
+    rme2 = BetaML.relative_mean_error(obs2,est2)
+    rmeb = BetaML.relative_mean_error(obsb,estb)
+    rme2b = BetaML.relative_mean_error(obs2b,est2b)
+
+    ms = BetaML.Scaler()
+    obss = BetaML.fit!(ms, obs)
+    ests = BetaML.predict(ms, est)
+    rmes =BetaML.relative_mean_error(obss, ests)
+
+    ms2 = BetaML.Scaler()
+    obs2s = BetaML.fit!(ms2, obs2)
+    est2s = BetaML.predict(ms2, est2)
+    rme2s = BetaML.relative_mean_error(obs2s, est2s)
+
+    msb = BetaML.Scaler()
+    obsbs = BetaML.fit!(msb, obsb)
+    estbs = BetaML.predict(msb, estb)
+    rmebs =BetaML.relative_mean_error(obsbs, estbs)
+
+    ms2 = BetaML.Scaler()
+    obs2s = BetaML.fit!(ms2, obs2)
+    est2s = BetaML.predict(ms2, est2)
+    rme2s = BetaML.relative_mean_error(obs2s, est2s)
+
+    ms2b = BetaML.Scaler()
+    obs2bs = BetaML.fit!(ms2b, obs2b)
+    est2bs = BetaML.predict(ms2b, est2b)
+    rme2bs = BetaML.relative_mean_error(obs2bs, est2bs)
+
+    masev_seq = BetaML.mase(obs,est;seq=true)
+    masev2_seq = BetaML.mase(obs2,est2;seq=true)
+    masevb_seq = BetaML.mase(obsb,estb;seq=true)
+    masev2b_seq = BetaML.mase(obs2b,est2b;seq=true)
+
+    @test masevb_seq > masev_seq
+
+    masev = BetaML.mase(obs,est;seq=false)
+    masev2 = BetaML.mase(obs2,est2;seq=false)
+    masevb = BetaML.mase(obsb,estb;seq=false)
+    masev2b = BetaML.mase(obs2b,est2b;seq=false)
+
+    @test masev ≈ masev2b ≈ rme2bs
+
+end
+
 # ==================================
 # New test
 println("** Testing pca()...")