Wrapping around Python PyStruct: structured SVM in Julia

external_model.py

from pystruct.models.base import StructuredModel
import numpy as np
import types

class ExternalModel(StructuredModel):
    """Interface definition for Structured Learners.
    This class defines what is necessary to use the structured svm.
    You have to implement at least joint_feature and inference.
    """
    def __init__(self, size_joint_feature=None, joint_feature=None,
                       inference=None, loss=None,
                       loss_augmented_inference=None, initialize=None):
        """Initialize the model.
        Needs to set self.size_joint_feature, the dimensionalty of the joint features for
        an instance with labeling (x, y).
        """
        self.size_joint_feature = size_joint_feature
        if joint_feature is not None:
            self.joint_feature = types.MethodType(joint_feature, self)
        if inference is not None:
            self.inference = types.MethodType(inference, self)
        if loss is not None:
            self.loss = types.MethodType(loss, self)
        if loss_augmented_inference is not None:
            self.loss_augmented_inference = types.MethodType(loss_augmented_inference, self)
        if initialize is not None:
            self.initialize = types.MethodType(initialize, self)

    def initialize(self, X, Y):
        pass

    def joint_feature(self, x, y):
        raise NotImplementedError()

    def inference(self, x, w, relaxed=None):
        raise NotImplementedError()

    def loss(self, y, y_hat):
        # hamming loss:
        return np.sum(y != y_hat)

    def loss_augmented_inference(self, x, y, w, relaxed=None):
        print("FALLBACK no loss augmented inference found")
        return self.inference(x, w)

regression_model.jl

# Reproducing Andrea Vedaldi's example on SSVM: regressing a real function
# example from http://www.robots.ox.ac.uk/~vedaldi/assets/svm-struct-matlab/tutorial/ssvm-tutorial.pdf
#
# Original description from Andrea Vedaldi:
#
#   Demonstrates using SSVM-STRUCT-MALTAB to learn a structured SVM
#   regressing a real function. The joint feature map is:
#
#     Psi(x,y) = [y yx yx^2 yx^3 -0.5y^2]'
#
#   which corresponds to learning a polynomial of the third order:
#
#     y(x;w) = argmax_y <Psi(x,y),w>
#            = w1/w5 + w2/w5 x + w4/w5 x^2 + w4/w5 x^3.
#
#   To avoid pathologies, y is constrained in the range [-1,1].
#
#   The loss is the absolute difference of the ground truth output
#   yi and the regressed output y(xi;w):
#
#     Delta(yi,y) = |y - yi|.
#
#   The marging rescaling formulation is used to construct the
#   surrogate loss:
#
#     L(w) = sup_y |y - yi| + <Psi(xi,y),w> - <Psi(xi,yi),w>


using PyCall; unshift!(PyVector(pyimport("sys")["path"]), "");

OneSlackSSVM = pyimport("pystruct.learners")["OneSlackSSVM"]  # learner
ExternalModel = pyimport("external_model")["ExternalModel"]   # model

println("Defining model functions...")
size_joint_feature = 5
function joint_feature(self, x, y)
    [y, y*x, y*x^2, y*x^3, -.5*y^2]
end

function inference(self, x, w; relaxed=None)
    z = (w[1] + w[2]*x + w[3]*x^2 + w[4]*x^3)
    if w[5] == 0.
        y_hat = sign(z)
    else
        y_all = [-1., 1., z/w[5]]
        y_hat = 0.
        loss = -Inf
        for c in y_all
            if -1. <= c <= 1.
                E = z*c - .5*c^2*w[5]
                if E > loss
                    loss = E
                    y_hat = c
                end
            end
        end
    end
    y_hat
end

function loss(self, y, y_hat)
    abs(y_hat - y)
end

function loss_augmented_inference(self, x, y, w; relaxed=None)
    z = (w[1] + w[2]*x + w[3]*x^2 + w[4]*x^3)
    if w[5] > 0
        y_all = [-1., 1., (z-1.)/w[5], (z+1.)/w[5]]
    else
        y_all = [-1., 1.]
    end
    y_hat = 0.
    augmentedloss = -Inf
    for c in y_all
        if -1. <= c <= 1.
            L = abs(c - y) + z*c - .5*c^2*w[5]
            if L > augmentedloss
                augmentedloss = L
                y_hat = c
            end
        end
    end
    y_hat
end

println("Creating model...")

juliamodel = pycall(ExternalModel, PyObject,
                    size_joint_feature,
                    joint_feature,
                    inference,
                    loss,
                    loss_augmented_inference)

x = collect(linspace(-pi,pi,21));
y = 0.5*sin(x);
y = y + 0.1*randn(size(y));

ssvm = pycall(OneSlackSSVM, PyObject, juliamodel)


println("Learning parameters...")
ssvm[:fit](x, y)

println("Evaluating results...")

w = ssvm[:w] # learned weights

minx, maxx = minimum(x), maximum(x)
xr = linspace(minx, maxx, 1024);
yr = linspace(-1,1,1024) ;

y_fit = ssvm[:predict](xr)

z = w[1] + w[2] * xr + w[3] * xr.^2 + w[4] * xr.^3;
F = yr*z' - 0.5 * yr.^2 * ones(size(z))' * w[5]; # scoring function

F_ = F .- maximum(F, 1); # column rescaled

println("Plotting results...")

using PyPlot

imshow(F_, extent=[minimum(x), maximum(x), -1, 1], origin="lower")
# plot(xr, y_fit)
plot(x, y, "o")
xlim([minimum(x), maximum(x)]);

Strange segmentation errors appear after the weight learning... Learning w should be enough for most use cases though. Strange.