MultivariateStats + StatsAPI

Project.toml

@@ -1,7 +1,7 @@
 name = "SpeciesDistributionToolkit"
 uuid = "72b53823-5c0b-4575-ad0e-8e97227ad13b"
 authors = ["Timothée Poisot <[email protected]>"]
-version = "0.0.10"
+version = "0.1.0"
 
 [deps]
 ArchGDAL = "c9ce4bd3-c3d5-55b8-8973-c0e20141b8c3"
@@ -17,6 +17,13 @@ SimpleSDMDatasets = "2c7d61d0-5c73-410d-85b2-d2e7fbbdcefa"
 SimpleSDMLayers = "2c645270-77db-11e9-22c3-0f302a89c64c"
 StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 
+[weakdeps]
+MultivariateStats = "6f286f6a-111f-5878-ab1e-185364afe411"
+StatsAPI = "82ae8749-77ed-4fe6-ae5f-f523153014b0"
+
+[extensions]
+MultivariateExtension = ["MultivariateStats", "StatsAPI"]
+
 [compat]
 ArchGDAL = "0.9, 0.10"
 Distances = "0.10"
@@ -30,4 +37,4 @@ Reexport = "1.2"
 SimpleSDMDatasets = "0.1"
 SimpleSDMLayers = "0.9"
 StatsBase = "0.33, 0.34"
-julia = "1.8"
+julia = "1.9"

ext/MultivariateExtension.jl

@@ -0,0 +1,63 @@
+module MultivariateExtension
+
+using SpeciesDistributionToolkit
+using MultivariateStats
+using StatsAPI
+
+# These types have a fit method
+types_to_fit = [:PCA, :PPCA, :KernelPCA, :Whitening, :MDS, :MetricMDS]
+
+# These types have a transform method
+types_to_transform = [:Whitening]
+
+# These types have a predict method
+types_to_predict = [:PCA, :PPCA, :KernelPCA, :MDS, :MetricMDS]
+
+for tf in types_to_fit
+    eval(
+        quote
+            function StatsAPI.fit(::Type{MultivariateStats.$tf}, X::Vector{T}; kwargs...) where {T <: SimpleSDMLayers.SimpleSDMLayer}
+                @assert SimpleSDMLayers._layers_are_compatible(X)
+                return StatsAPI.fit(MultivariateStats.$tf, Array(X); kwargs...)
+            end
+        end
+    )
+end
+
+for tf in types_to_transform
+    eval(
+        quote
+            function MultivariateStats.transform(f::MultivariateStats.$tf, X::Vector{T}; kwargs...) where {T <: SimpleSDMLayers.SimpleSDMLayer}
+                @assert SimpleSDMLayers._layers_are_compatible(X)
+                D = MultivariateStats.transform(f, Array(X); kwargs...)
+                O = [similar(X[1], eltype(D)) for i in 1:length(X)]
+                for i in axes(O, 1)
+                    for (j, k) in enumerate(keys(O[i]))
+                        O[i][k] = D[i, j]
+                    end
+                end
+                return O
+            end
+        end
+    )
+end
+
+for tf in types_to_predict
+    eval(
+        quote
+            function StatsAPI.predict(f::MultivariateStats.$tf, X::Vector{T}; kwargs...) where {T <: SimpleSDMLayers.SimpleSDMLayer}
+                @assert SimpleSDMLayers._layers_are_compatible(X)
+                D = StatsAPI.predict(f, Array(X); kwargs...)
+                O = [similar(X[1], eltype(D)) for i in 1:MultivariateStats.outdim(M)]
+                for i in axes(O, 1)
+                    for (j, k) in enumerate(keys(O[i]))
+                        O[i][k] = D[i, j]
+                    end
+                end
+                return O
+            end
+        end
+    )
+end
+
+end

src/SpeciesDistributionToolkit.jl

@@ -22,6 +22,9 @@ using Reexport
 @reexport using Fauxcurrences
 @reexport using Phylopic
 
+# Quality of life functions
+include("quality_of_life.jl")
+
 # SimpleSDMLayers to wrap everything together
 include("integrations/datasets_layers.jl")
 
@@ -31,7 +34,7 @@ include("integrations/gbif_layers.jl")
 # GBIF and Phylopic integration
 include("integrations/gbif_phylopic.jl")
 
-# Plotting
+# Makie recipes
 include("integrations/makie.jl")
 
 # Functions for IO

src/integrations/makie.jl

@@ -1,4 +1,3 @@
-# Function to turn a layer into something (Geo)Makie can use
 function sprinkle(layer::T) where {T <: SimpleSDMLayer}
     return (
         longitudes(layer),

src/quality_of_life.jl

@@ -0,0 +1,12 @@
+function Base.Array(layers::Vector{T}) where {T <: SimpleSDMLayers.SimpleSDMLayer}
+    @assert SimpleSDMLayers._layers_are_compatible(layers)
+    k = reduce(intersect, [findall(!isnothing, grid(layer)) for layer in layers])
+    R = SimpleSDMLayers._inner_type(first(layers))
+    X = Matrix{R}(undef, length(layers), length(k))
+    for i in eachindex(k)
+        for j in eachindex(layers)
+            X[j,i] = layers[j][k[i]]
+        end
+    end
+    return X
+end