Commit 48a51446 authored by Corson N. Areshenkoff's avatar Corson N. Areshenkoff

Add parallel transport

parent 36a0dbf0
......@@ -10,6 +10,7 @@ export(spd.expmap)
export(spddot)
export(spd.whiten)
export(spd.correlation)
export(spd.transport)
importFrom(expm,expm)
importFrom(expm,logm)
......
logm2 <- function(x){
tryCatch(logm(x, method = 'Eigen'),
error = function(e) logm(x))
lx <- logm(x)
return((x + t(x))/2)
}
expm2 <- function(x){
tryCatch(expm(x, method = 'R_Eigen'),
error = function(e) expm(x))
ex <- expm(x)
return((x + t(x))/2)
}
#' Projection onto the tangent space
#' Projection from the tangent space
#'
#' Function projects a tangent vector x (a symmetric matrix) at a point p
#' onto the space of SPD matrix.
#' onto the space of SPD matrices.
#'
#' @param x A symmetric matrix
#' @param p The point to whose tangent space x belongs
......
......@@ -5,7 +5,7 @@
#' \code{t=1} returns \code{y}, and \code{t=.5} returns the mean of \code{x} and \code{y}
#'
#' @param x,y Symmetric, positive-definite matrices
#' @param t Interpolation parameter in \code{[0,1]}
#' @param t Interpolation parameter in \code{[0,infinity)}
#' @param method Type of interpolation (see details)
#' @details Allowable distance measures are
#' \itemize{
......@@ -19,8 +19,8 @@ spd.interpolate <- function(x, y, t, method = 'euclidean', ...){
if (!'spd.mat' %in% input.type(x) | !'spd.mat' %in% input.type(y)){
stop('Both inputs must be positive definite matrices')
}
if (t < 0 | t > 1){
stop('t must lie in the interval [0,1]')
if (t < 0){
stop('t must be greater than 1')
}
# Wrapper
......
......@@ -15,7 +15,7 @@
#' \item{"logeuclidean": }{Computed by taking the arithmetic mean of the logarithms
#' of the matrices in \code{x}, and then projecting back onto the space of spd matrices.
#' In general, better behaved than the arithmetic mean.}
#' \item{"riemannian": }{The Riemmanian p-mean by. Smoothly interpolates between the
#' \item{"riemannian": }{The Riemmanian p-mean. Smoothly interpolates between the
#' harmonic mean at \code{p = -1} to the geometric mean at \code{p = 0} to the
#' arithmetic mean at \code{p = 1}.
#' Is approximated iteratively using the fixed point algorithm described by
......
#' Parallel transport of a tangent vector
#'
#' Whitens an SPD matrix using the procedure advocated by Ng, et al. (2014)
#' Parallel transports a vector \code{x} in the tangent space at \code{"from"} to the
#' tangent space at \code{"to"} using the Schild's ladder algorithm.
#'
#' @param x An SPD matrix to be whitened
#' @param p The SPD baseline whitening matrix
#' @param unwhiten Logical. If TRUE, reverses a previously applied whitening transform.
#' @return A symmetric, positive-definite matrix.
to <- cov(matrix(rnorm(30), ncol = 3))
from <- cov(matrix(rnorm(30), ncol = 3))
x <- cov(matrix(rnorm(30), ncol = 3))
x <- spd.logmap(x, from)
#' @param x A tangent vector (symmetric matrix)
#' @param to The SPD matrix to whose tangent space to move x
#' @param from The SPD matrix to whose tangent space x belongs
#' @param nsteps The number of steps in the geodesic connecting to and from
#' @return A symmetric matrix -- a tangent vector at to.
spd.transport <- function(x, to, from, nsteps = 10){
......@@ -25,11 +22,19 @@ spd.transport <- function(x, to, from, nsteps = 10){
}
# Generate discrete geodesic from source to target
geod <- lapply(seq(0, 1, length.out = nsteps), function(i) {
geod <- lapply(seq(0, 1, length.out = nsteps+1)[-1], function(i) {
spd.interpolate(from, to, t = i, method = 'riemannian')
})
return(NULL)
# Iterate over rungs
trans <- spd.expmap(x, from)
source <- from
for (i in 1:nsteps){
m <- spd.interpolate(trans, geod[[i]], t = .5, method = 'riemannian')
trans <- spd.interpolate(source, m, t = 2, method = 'riemannian')
source <- geod[[i]]
}
return(spd.logmap(trans, p = to))
}
#' Vectorize and unvectorize a matrix
#'
#' Function implements several distance measures between symmetric,
#' positive-definite matrices, not all of which are proper metrics.
#' Function converts between a square symmetric matrix and a vector of the
#' lower triangular elements + diagonal. Allows optional scaling of the off-diagonal
#' entries in order to preserve the norm.
#'
#' @param x Either a square matrix, or a vector contains the lower triangular
#' elements of x (including the diagonal elements)
......
......@@ -2,7 +2,7 @@
% Please edit documentation in R/spd-expmap.R
\name{spd.expmap}
\alias{spd.expmap}
\title{Projection onto the tangent space}
\title{Projection from the tangent space}
\usage{
spd.expmap(x, p)
}
......@@ -16,5 +16,5 @@ A symmetric, positive-definite matrix.
}
\description{
Function projects a tangent vector x (a symmetric matrix) at a point p
onto the space of SPD matrix.
onto the space of SPD matrices.
}
......@@ -9,7 +9,7 @@ spd.interpolate(x, y, t, method = "euclidean", ...)
\arguments{
\item{x, y}{Symmetric, positive-definite matrices}
\item{t}{Interpolation parameter in \code{[0,1]}}
\item{t}{Interpolation parameter in \code{[0,infinity)}}
\item{method}{Type of interpolation (see details)}
}
......
......@@ -4,7 +4,7 @@
\alias{spd.logmap}
\title{Projection onto the tangent space}
\usage{
spd.logmap(x, p)
spd.logmap(x, p = NULL)
}
\arguments{
\item{x}{A symmetric positive definite matrix}
......
......@@ -30,7 +30,7 @@ matrices. Several methods are implemented:
\item{"logeuclidean": }{Computed by taking the arithmetic mean of the logarithms
of the matrices in \code{x}, and then projecting back onto the space of spd matrices.
In general, better behaved than the arithmetic mean.}
\item{"riemannian": }{The Riemmanian p-mean by. Smoothly interpolates between the
\item{"riemannian": }{The Riemmanian p-mean. Smoothly interpolates between the
harmonic mean at \code{p = -1} to the geometric mean at \code{p = 0} to the
arithmetic mean at \code{p = 1}.
Is approximated iteratively using the fixed point algorithm described by
......
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/spd-transport.R
\name{spd.transport}
\alias{spd.transport}
\title{Parallel transport of a tangent vector}
\usage{
spd.transport(x, to, from, nsteps = 10)
}
\arguments{
\item{x}{A tangent vector (symmetric matrix)}
\item{to}{The SPD matrix to whose tangent space to move x}
\item{from}{The SPD matrix to whose tangent space x belongs}
\item{nsteps}{The number of steps in the geodesic connecting to and from}
}
\value{
A symmetric matrix -- a tangent vector at to.
}
\description{
Parallel transports a vector \code{x} in the tangent space at \code{"from"} to the
tangent space at \code{"to"} using the Schild's ladder algorithm.
}
......@@ -14,8 +14,9 @@ elements of x (including the diagonal elements)}
the norm of the resulting vector}
}
\description{
Function implements several distance measures between symmetric,
positive-definite matrices, not all of which are proper metrics.
Function converts between a square symmetric matrix and a vector of the
lower triangular elements + diagonal. Allows optional scaling of the off-diagonal
entries in order to preserve the norm.
}
\details{
If input is a square matrix, converts the lower triangular elements
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment