Add orthogonal and semi-orthogonal matrix types

Abierto
  • ¿Algo no está bien?
#3,279 3 comentarios 0 reacciones 0 asignados Ver en GitHub

Los mantenedores suelen responder en 1 día

Nadie ha tomado este issue todavía.

Evaluación

Dificultad
5/5
Tiempo estimado
Más de una semana
Aptitud para principiantes
35/100
Tipo de issue
Nueva funcionalidad
Claridad
Bastante claro
Estado de actividad
Estancado
Stack tecnológico
cpp
Área
backend

Línea de trabajo

Empieza leyendo las discusiones enlazadas del foro de Stan, la tesis de Michael Jauch y la implementación de stan_semiorthogonal_transforms. Compara los tipos propuestos orthogonal_matrix[K] y semi_orthogonal_matrix[K, N] con la parametrización mediante reflectores y sus reglas de conteo de parámetros. Se considera terminado cuando el proyecto cuenta con una implementación acordada que cubre los tipos solicitados y operaciones eficientes de multiplicación, determinante y resolución de sistemas lineales.

Escrito por el modelo de indexación a partir del texto del issue.

Descripción

Add orthogonal and semi-orthogonal matrix types. These are described in https://en.wikipedia.org/wiki/Orthogonal_matrix. There are numerous applications and many of them are detailed in Michael Jauch's PhD thesis https://dukespace.lib.duke.edu/items/da7e8aed-efe9-4006-bbd9-ecc3136e0d3b.

An efficient parameterization is given in https://github.com/sethaxen/stan_semiorthogonal_transforms and the main Stan function file https://github.com/sethaxen/stan_semiorthogonal_transforms/blob/main/transforms/semiorthogonal_reflector_transform.stan. These are compared to other paramerterizations in https://discourse.mc-stan.org/t/parameterizing-an-orthonormal-matrix-the-stiefel-manifold/33211/11?u=spinkney.

The proposal is to have two types

parameters {
 orthogonal_matrix[K] x;
 semi_orthogonal_matrix[K, N] y;
}

The orthogonal matrix is a square K x K matrix with orthonormal

The number of parameters is in the above links as

/**
 * Return the number of unconstrained parameters needed for a N x K semiorthogonal matrix Q.
 * 
 * @param N number of rows of Q
 * @param K number of columns of Q (must be <= N)
 * @param special whether Q is special orthogonal (i.e. det(Q)=+1). Must be in (0,1). Only
 *                checked if N==K.
 * @return nparams Number of unconstrained parameters
 */
int semiorthogonal_reflector_num_params(int N, int K, int special){
  return N * K - (K * (K - 1)) %/% 2 - (N==K && special);
}

This contrasts with the "naive" method of generating semi-orthogonal matrices of K x N. The Householder reflector parameterization has choose(K, 2) fewer parameters. In the case of generating a square orthogonal matrix there is an additional reduction in the number of parameters by 1 due to restricting diagonal sign to be positive. We don't need the "special" keyword as we would fix the sign indeterminancy to positive 1, just like we currently do with the QR decomposition in Stan.

As noted in the forum post by Seth that the implementation would be adding a type and then having efficient ways to do multiplication, determinant, and linear solves with these factors:

Finally, reflectors_factors imagines what if Stan added a native semiorthogonal_matrix[K,N] data type that stored the reflector factors and for which certain operations like multiplication, determinant, and linear solve were implemented in terms of those factors. Efficient algorithms exist for each of these and are implemented in Eigen and LAPACK. This is really the only option with these Householder-inspired transforms if one needs to efficiently infer a large orthogonal matrix.

Adding various questions and applications about orthogonal matrices from the forums

Lenguaje dominante
C++
Estrellas
845
Forks
225
Merge medio
3 d 9 h
PR fusionados (30 d)
12

Preparar el entorno

Primeros pasos

  1. Lee el issue completo y luego la guía de contribución del proyecto.
  2. Comenta en el issue que vas a ocuparte — evita que dos personas hagan lo mismo.
  3. Haz un fork del repositorio y trabaja en una rama.
  4. Abre un pull request que haga referencia al número del issue.

Más de stan-dev/math

Todos los issues de stan-dev/math

Issues similares

Más issues de C++

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.