HDBSCAN and DBCV example test application in Julia

hdbscan_example_hyperparameter_optimization.jl

# ==========================================================
# HDBSCAN hyperparameter optimization using DBCV validation
# ==========================================================

using Pkg

# ----------------------------------------------------------
# Install missing packages
# ----------------------------------------------------------

function ensure_package(pkg::String)
    try
        Base.eval(Main, Meta.parse("using $pkg"))
        println("$pkg already available")
    catch
        println("Installing $pkg ...")
        Pkg.add(pkg)
    end
end

for pkg in [
    "Distributions",
    "HDBSCAN",
    "DensityBasedClusteringValidation",
    "Plots"
]
    ensure_package(pkg)
end


# ----------------------------------------------------------
# Load packages
# ----------------------------------------------------------

using Random
using LinearAlgebra
using Distributions
using HDBSCAN
using DensityBasedClusteringValidation
using Plots
using Printf


# ----------------------------------------------------------
# Generate synthetic data
# ----------------------------------------------------------

Random.seed!(42)

n = 500
omega = 0.25
Sigma = omega^2 * I

d1 = MvNormal([-1.0, 1.0], Sigma)
d2 = MvNormal([-0.5, 0.5], Sigma)

# 1000 samples × 2 features
X = vcat(
    rand(d1, n)',
    rand(d2, n)'
)

println("Data size: ", size(X))


# ----------------------------------------------------------
# HDBSCAN optimization function
# ----------------------------------------------------------

function optimize_hdbscan(X)

    # Hyperparameter grid
    min_cluster_sizes = [3, 4, 5, 6, 10, 15, 20, 30, 50, 60, 70, 80, 90, 100]
    min_samples_values = [5, 10, 12, 15, 17, 20, 30, 35, 40]

    best_dbcv = -Inf
    best_min_cluster_size = nothing
    best_min_samples = nothing
    best_labels = nothing
    best_model = nothing


    for min_cluster_size in min_cluster_sizes

        for min_samples in min_samples_values

            println("\n--------------------------------")
            println(
                "Testing min_cluster_size = ",
                min_cluster_size,
                ", min_samples = ",
                min_samples
            )


            # -----------------------------
            # HDBSCAN
            # -----------------------------

            model = Hdbscan(
                min_cluster_size,
                min_samples;
                metric = "euclidean",
                copy = true
            )

            fit!(model, X)

            labels = HDBSCAN.labels(model)


            # Skip useless solutions
            if all(labels .== -1)
                println("All points classified as noise")
                continue
            end


            # -----------------------------
            # DBCV score
            # -----------------------------

            dbcv_score = dbcv(
                X,
                labels
            )


            println(
                "Clusters = ",
                HDBSCAN.nclusters(model),
                ", DBCV = ",
                @sprintf("%.10f", dbcv_score)
            )


            # Save best model

            if dbcv_score > best_dbcv

                best_dbcv = dbcv_score
                best_min_cluster_size = min_cluster_size
                best_min_samples = min_samples
                best_labels = labels
                best_model = model

            end

        end
    end


    return (
        best_dbcv,
        best_min_cluster_size,
        best_min_samples,
        best_labels,
        best_model
    )

end


# ----------------------------------------------------------
# Run optimization
# ----------------------------------------------------------

(
    best_dbcv,
    best_min_cluster_size,
    best_min_samples,
    best_labels,
    best_model
) = optimize_hdbscan(X)



# ----------------------------------------------------------
# Print best result
# ----------------------------------------------------------

println("\n================================")
println("BEST HDBSCAN CONFIGURATION")
println("================================")

println(
    "min_cluster_size = ",
    best_min_cluster_size
)

println(
    "min_samples = ",
    best_min_samples
)

println(
    "Number of clusters = ",
    HDBSCAN.nclusters(best_model)
)

println(
    "DBCV score = ",
    @sprintf("%.10f", best_dbcv)
)


#=

# ----------------------------------------------------------
# Plot best clustering
# ----------------------------------------------------------

scatter(
    X[:,1],
    X[:,2],
    marker_z = best_labels,
    markersize = 4,
    legend = false,
    xlabel = "Feature 1",
    ylabel = "Feature 2",
    title = "Best HDBSCAN clustering (DBCV=$(round(best_dbcv,digits=5)))"
)=#