Merge pull request #16 from tiemvanderdeure/statsapi_predict

use StatsAPI.predict, add mlj metadata, better tests, some other minor changes

Author: tiemvanderdeure

docs/src/usage/quickstart.md

@@ -26,7 +26,7 @@ bradypus_model = maxnet(p_a, env)
 prediction = predict(bradypus_model, env)
 ```
 
-There are numerous settings that can be tweaked to change the model fit. These are documentated in the documentatoin for the `maxnet`(@ref) and `Maxnet.predict`(@ref) functions.
+There are numerous settings that can be tweaked to change the model fit. These are documentated in the documentation for the `maxnet`(@ref) and `predict`(@ref) functions.
 
 ### Model settings
 The two most important settings to change when running Maxnet is the feature classes selected and the regularization factor.

src/Maxnet.jl

@@ -4,12 +4,13 @@ import Tables, Lasso, GLMNet, Interpolations, CategoricalArrays, GLM, SparseArra
 import StatsAPI, StatsBase, Statistics
 import MLJModelInterface as MMI
 
+using StatsAPI: predict
 using GLM: IdentityLink, CloglogLink, LogitLink, LogLink
 using MLJModelInterface: Continuous, Binary, Multiclass, Count
 
 export IdentityLink, CloglogLink, LogitLink, LogLink # re-export relevant links
 export maxnet, predict, complexity
-export LinearFeature, CategoricalFeature, QuadraticFeature, ProductFeature, ThresholdFeature, HingeFeature
+export LinearFeature, CategoricalFeature, QuadraticFeature, ProductFeature, ThresholdFeature, HingeFeature, AbstractFeatureClass
 export MaxnetBinaryClassifier
 
 

src/maxnet_function.jl

@@ -43,7 +43,7 @@ Variables selected: [:frs6190_ann, :h_dem, :pre6190_l1, :pre6190_l10, :pre6190_l
 function maxnet(
     presences::BitVector, predictors; 
     features = default_features(sum(presences)),
-    regularization_multiplier::Float64 = 1.0,
+    regularization_multiplier = 1.0,
     regularization_function = default_regularization,
     addsamplestobackground::Bool = true, weight_factor::Float64 = 100.,
     n_knots::Int = 50,
@@ -82,7 +82,7 @@ function _maxnet(
     presences::BitVector, 
     predictors, 
     features::Vector{<:AbstractFeatureClass},
-    regularization_multiplier::Float64,
+    regularization_multiplier,
     regularization_function,
     addsamplestobackground::Bool, 
     weight_factor::Float64,

src/mlj_interface.jl

@@ -46,12 +46,24 @@ end
 """
 MaxnetBinaryClassifier
 
-MMI.input_scitype(::Type{<:MaxnetBinaryClassifier}) =
-        MMI.Table{<:Union{<:AbstractVector{<:Continuous}, <:AbstractVector{<:Multiclass}}} #{<:Union{<:Continuous <:Multiclass}}
-    
-MMI.target_scitype(::Type{<:MaxnetBinaryClassifier}) = AbstractVector{Multiclass{2}}# AbstractVector{<:MMI.Finite}
+MMI.metadata_pkg(
+    MaxnetBinaryClassifier;
+    name = "Maxnet",
+    uuid = "81f79f80-22f2-4e41-ab86-00c11cf0f26f",
+    url = "https://github.com/tiemvanderdeure/Maxnet.jl",
+    is_pure_julia = false,
+    package_license = "MIT",
+    is_wrapper = false    
+)
 
-MMI.load_path(::Type{<:MaxnetBinaryClassifier}) = "Maxnet.MaxnetBinaryClassifier"
+MMI.metadata_model(
+    MaxnetBinaryClassifier;
+    input_scitype = MMI.Table(MMI.Continuous, MMI.Finite),
+    target_scitype = AbstractVector{<:MMI.Finite{2}},
+    load_path = "Maxnet.MaxnetBinaryClassifier",
+    human_name = "Maxnet",
+    reports_feature_importances=false
+)
 
 function MMI.fit(m::MaxnetBinaryClassifier, verbosity::Int, X, y)
     # convert categorical to boolean

src/predict.jl

@@ -22,10 +22,10 @@ A `Vector` with the resulting predictions.
 using Maxnet
 p_a, env = Maxnet.bradypus();
 bradypus_model = maxnet(p_a, env; features = "lq")
-prediction = Maxnet.predict(bradypus_model, env)
+prediction = predict(bradypus_model, env)
 ```
 """
-function predict(m::MaxnetModel, x; link = CloglogLink(), clamp = false)
+function StatsAPI.predict(m::MaxnetModel, x; link = CloglogLink(), clamp = false)
     predictors = Tables.columntable(x)
     for k in keys(m.predictor_data)
         k in keys(predictors) || error("$k is not found in the predictors")

src/types.jl

@@ -1,5 +1,5 @@
 struct MaxnetModel 
-    path::Union{GLMNet.GLMNetPath, Lasso.LassoPath}
+    path::GLMNet.GLMNetPath
     features::Vector{<:AbstractFeatureClass}
     columns::Vector{ModelMatrixColumn}
     coefs::AbstractVector
@@ -11,7 +11,7 @@ struct MaxnetModel
 end
 
 function Base.show(io::IO, mime::MIME"text/plain", m::MaxnetModel)
-    vars_selected = mapreduce(Maxnet._var_keys, (x, y) -> unique(vcat(x, y)), selected_features(m))
+    vars_selected = mapreduce(Maxnet._var_keys, (x, y) -> unique(vcat(x, y)), selected_features(m); init = Symbol[])
 
     println(io, "Fit Maxnet model")
 

test/runtests.jl

@@ -1,4 +1,4 @@
-using Maxnet, Test, MLJBase, Statistics
+using Maxnet, Test, Statistics
 
 p_a, env = Maxnet.bradypus()
 env1 = map(e -> [e[1]], env) # just the first row
@@ -14,42 +14,81 @@ env1 = map(e -> [e[1]], env) # just the first row
     @test Maxnet.features_from_string("lqpt") == [LinearFeature(), CategoricalFeature(), QuadraticFeature(), ProductFeature(), ThresholdFeature()]
 
     @test Maxnet.default_features(100) == [LinearFeature(), CategoricalFeature(), QuadraticFeature(), HingeFeature(), ProductFeature()]
+
+    @test Maxnet.hinge(1:5, 3) ==    [ 
+     #  1:5   3:5  1:3  1:5 
+        0.0   0.0  0.0  0.0
+        0.25  0.0  0.5  0.25
+        0.5   0.0  1.0  0.5
+        0.75  0.5  1.0  0.75
+        1.0   1.0  1.0  1.0
+    ]
+    @test size(Maxnet.hinge(1:200)) == (200, 98)
+
+    presence, predictors = Maxnet.addsamples([true, false, false], (a = [1,2,3], b = [1,2,3]))
+    @test presence == [true, false, false, false]
+    @test predictors == (a = [1,2,3,1], b = [1,2,3,1])
 end
 
 @testset "Maxnet" begin
     # some class combinations and keywords
-    m = Maxnet.maxnet(p_a, env; features = "lq");
-    Maxnet.maxnet(p_a, env; features = "lqp", regularization_multiplier = 2.);
-    Maxnet.maxnet(p_a, env; features = "lqh", regularization_multiplier = 5., nknots = 10);
-    Maxnet.maxnet(p_a, env; features = "lqph", weight_factor = 10.);
+    m = maxnet(p_a, env; features = "lq");
+    m2 = maxnet(p_a, env)
+    m3 = maxnet(p_a, env[(:cld6190_ann, :h_dem)])
+    m4 = maxnet(p_a, env[(:ecoreg,)], addsamplestobackground =false)
+    m5 = maxnet(p_a, env[(:cld6190_ann, :h_dem)]; features = "ht", n_knots = 3)
 
-    # test the result
+    # test the results
     @test m.entropy ≈ 6.114650341746531
     @test complexity(m) == 21
-
+    @test m2.features == [LinearFeature(), CategoricalFeature(), QuadraticFeature(), HingeFeature(), ProductFeature()]
+    @test m3.features == [LinearFeature(), QuadraticFeature(), HingeFeature(), ProductFeature()]
+    @test m4.features == [CategoricalFeature()]
+    @test m5.features == [HingeFeature(), ThresholdFeature()]
+    @test length(m5.columns) == 14 # (n-1)*2 hinge columns and n threshold columns for each variable
+    
     # predictions
-    prediction = Maxnet.predict(m, env)
+    prediction = predict(m, env)
     @test Statistics.mean(prediction[p_a]) > Statistics.mean(prediction[.~p_a])
     @test minimum(prediction) > 0.
     @test maximum(prediction) < 1.
     @test mean(prediction) ≈ 0.24375837576014572 atol=1e-4
 
     # check that clamping works
     # clamp shouldn't change anything in this case
-    @test prediction == Maxnet.predict(m, env; clamp = true)
+    @test prediction == predict(m, env; clamp = true)
 
     # predict with a crazy extrapolation
     env1_extrapolated = merge(env1, (;cld6190_ann = [100_000]))
     env1_max_cld = merge(env1, (;cld6190_ann = [maximum(env.cld6190_ann)]))
 
     # using clamp the prediction uses the highest cloud
-    @test Maxnet.predict(m, env1_extrapolated; link = IdentityLink(), clamp = true) == 
-        Maxnet.predict(m, env1_max_cld; link = IdentityLink()) 
+    @test predict(m, env1_extrapolated; link = IdentityLink(), clamp = true) == 
+        predict(m, env1_max_cld; link = IdentityLink()) 
+
+    # test that maxnet works if no features are selected
+    empty_model = maxnet(p_a, env; regularization_multiplier = 1000);
+    @test complexity(empty_model) == 0
+    @test Maxnet.selected_features(empty_model) == Symbol[]
+    @test length(unique(predict(empty_model, env))) == 1
 end
 
 @testset "MLJ" begin
+    using MLJBase
     mn = Maxnet.MaxnetBinaryClassifier
 
+    # Test model metadata
+    @test name(mn) == "MaxnetBinaryClassifier"
+    @test human_name(mn) == "Maxnet"
+    @test package_name(mn) == "Maxnet"
+    @test !supports_weights(mn)
+    @test !is_pure_julia(mn)
+    @test is_supervised(mn)
+    @test package_license(mn) == "MIT"
+    @test prediction_type(mn) == :probabilistic
+    @test input_scitype(mn) == Table{<:Union{AbstractVector{<:Continuous}, AbstractVector{<:Finite}}}
+    @test hyperparameters(mn) == (:features, :regularization_multiplier, :regularization_function, :weight_factor, :link, :clamp, :kw)
+
     # convert to continuous
     cont_keys = collect(key => Continuous for key in keys(env) if key !== :ecoreg)
     env_typed = MLJBase.coerce(env, cont_keys...)
@@ -64,14 +103,13 @@ end
     # make the equivalent model without mlj
     model = Maxnet.maxnet((p_a), env_typed; features = "lqph");
 
-
     # predict via MLJBase
     mljprediction = MLJBase.predict(mach2, env_typed)
     mlj_true_probability = pdf.(mljprediction, true)
 
     # test that this predicts the same as the equivalent model without mlj
 
-    @test all(Maxnet.predict(model, env_typed) .≈ mlj_true_probability)
+    @test all(predict(model, env_typed) .≈ mlj_true_probability)
 
     @test Statistics.mean(mlj_true_probability[p_a]) > Statistics.mean(mlj_true_probability[.~p_a])
     @test minimum(mlj_true_probability) > 0.