Prototype for alternating gradient descent with Optim.jl

altgradpca.jl

using Optim

# low-dimensional data embedded in high-dimension
data = randn(100,5)*randn(5,100);

# a container for the parameters we fit
immutable PCA{T<:Real}
    X::Matrix{T}
    Y::Matrix{T}
end
PCA(m::Int,n::Int,r::Int) = PCA(randn(m,r),randn(r,n))

# initialize model
const p = PCA(100,100,5)

# total objective function is:
#       f(X,Y) = f1(X) + f2(Y)
f1(X) = 0.5*sumabs2(data - X*p.Y)
f2(Y) = 0.5*sumabs2(data - p.X*Y)

# We will alternate between updating X and Y since 
# the subproblems are convex and this may lead to
# better/quicker convergence

d1 = DifferentiableFunction(f1)
d2 = DifferentiableFunction(f2)

gd = GradientDescent()
o = OptimizationOptions()

s1 = Optim.initial_state(gd, o, d1, p.X)
s2 = Optim.initial_state(gd, o, d2, p.Y)

for iter = 1:iterations
    Optim.update!(d1, s1, gd)
    Optim.update!(d2, s2, gd)
end