K-Means Clustering (MSSC)¶

using System;
using System.IO;
using System.Globalization;
using Hexaly.Optimizer;

public class Kmeans : IDisposable
{
    // Data properties
    int nbObservations;
    int nbDimensions;
    int k;

    double[][] coordinatesData;

    // Hexaly Optimizer
    HexalyOptimizer optimizer;

    // Decisions
    HxExpression[] clusters;

    // Objective
    HxExpression obj;

    public Kmeans(int k)
    {
        optimizer = new HexalyOptimizer();
        this.k = k;
    }

    // Read instance data
    public void ReadInstance(string fileName)
    {
        using (StreamReader input = new StreamReader(fileName))
        {
            string[] splittedLine = input.ReadLine().Split();

            nbObservations = int.Parse(splittedLine[0]);
            nbDimensions = int.Parse(splittedLine[1]);

            coordinatesData = new double[nbObservations][];
            for (int o = 0; o < nbObservations; ++o)
            {
                splittedLine = input.ReadLine().Split();
                coordinatesData[o] = new double[nbDimensions];
                for (int d = 0; d < nbDimensions; ++d)
                    coordinatesData[o][d] = double.Parse(
                        splittedLine[d],
                        CultureInfo.InvariantCulture
                    );
            }
        }
    }

    public void Dispose()
    {
        if (optimizer != null)
            optimizer.Dispose();
    }

    public void Solve(int limit)
    {
        // Declare the optimization model
        HxModel model = optimizer.GetModel();

        // Set decisions: clusters[c] represents the points in cluster c
        clusters = new HxExpression[k];
        for (int c = 0; c < k; ++c)
            clusters[c] = model.Set(nbObservations);

        // Each point must be in one cluster and one cluster only
        model.Constraint(model.Partition(clusters));

        // Coordinates of points
        HxExpression coordinates = model.Array(coordinatesData);

        // Compute variances
        HxExpression[] variances = new HxExpression[k];
        for (int c = 0; c < k; ++c)
        {
            HxExpression cluster = clusters[c];
            HxExpression size = model.Count(cluster);

            // Compute the centroid of the cluster
            HxExpression centroid = model.Array();
            for (int d = 0; d < nbDimensions; ++d)
            {
                HxExpression coordinateLambda = model.LambdaFunction(
                    o => model.At(coordinates, o, model.CreateConstant(d))
                );
                centroid.AddOperand(
                    model.If(size == 0, 0, model.Sum(cluster, coordinateLambda) / size)
                );
            }

            // Compute the variance of the cluster
            HxExpression variance = model.Sum();
            for (int d = 0; d < nbDimensions; ++d)
            {
                HxExpression dimensionVarianceLambda = model.LambdaFunction(
                    o =>
                        model.Pow(
                            model.At(coordinates, o, model.CreateConstant(d))
                                - model.At(centroid, model.CreateConstant(d)),
                            2
                        )
                );
                HxExpression dimensionVariance = model.Sum(cluster, dimensionVarianceLambda);
                variance.AddOperand(dimensionVariance);
            }
            variances[c] = variance;
        }

        // Minimize the total variance
        obj = model.Sum(variances);
        model.Minimize(obj);

        model.Close();

        // Parametrize the optimizer
        optimizer.GetParam().SetTimeLimit(limit);

        optimizer.Solve();
    }

    /* Write the solution in a file in the following format:
     *  - objective value
     *  - k
     *  - for each cluster, a line with the elements in the cluster (separated by spaces) */
    public void WriteSolution(string fileName)
    {
        using (StreamWriter output = new StreamWriter(fileName))
        {
            output.WriteLine(obj.GetDoubleValue());
            output.WriteLine(k);
            for (int c = 0; c < k; ++c)
            {
                HxCollection clusterCollection = clusters[c].GetCollectionValue();
                for (int i = 0; i < clusterCollection.Count(); ++i)
                    output.Write(clusterCollection[i] + " ");
                output.WriteLine();
            }
            output.Close();
        }
    }

    public static void Main(string[] args)
    {
        if (args.Length < 1)
        {
            Console.WriteLine("Usage: Kmeans inputFile [outputFile] [timeLimit] [k value]");
            Environment.Exit(1);
        }

        string instanceFile = args[0];
        string outputFile = args.Length > 1 ? args[1] : null;
        string strTimeLimit = args.Length > 2 ? args[2] : "5";
        string k = args.Length > 3 ? args[3] : "2";
        using (Kmeans model = new Kmeans(int.Parse(k)))
        {
            model.ReadInstance(instanceFile);
            model.Solve(int.Parse(strTimeLimit));
            if (outputFile != null)
                model.WriteSolution(outputFile);
        }
    }
}