Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
projected move regression
using LinearAlgebra
function linreg(xs, ys)
# Axs + b ≈ ys
@assert size(xs,2) == size(ys,2)
nobs = size(xs, 2)
xs_ = [xs; transpose(ones(nobs))]
ys_ = [ys; transpose(ones(nobs))]
A_ = ys_ / xs_
@assert A_[end, end] ≈ 1
@assert norm(A_[end, 1:end-1]) < sqrt(eps(eltype(A_)))
A = A_[1:end-1, 1:end-1]
b = A_[1:end-1, end]
A, b
end
function projected_move_regression(proj, xs, ys)
# ys ≈ proj * R * xs .+ b
# where R is a rotation
projR, b = linreg(xs, ys)
A = proj \ projR
U,S,V = svd(A)
Vt = transpose(V)
@assert A ≈ U * Diagonal(S) * Vt
R = U * Vt
R, b
end
using Rotations, StaticArrays
using Test
A = randn(2,2)
b = randn(2)
xs = randn(2,10)
ys = A * xs .+ b
A_reg, b_reg = linreg(xs, ys)
@test A_reg ≈ A
@test b_reg ≈ b
rot_truth = rand(RotMatrix{3})
proj_trans_truth = @SVector randn(2)
proj = [1 0 0; 0 1 0]
N = 5
xs = randn(3,N)
ys = proj * rot_truth * xs .+ proj_trans_truth
R, b = projected_move_regression(proj, xs, ys)
@test R ≈ rot_truth
@test b ≈ proj_trans_truth
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.