spnorm: Spherical Normal Distribution In Riemann: Learning with Data on Riemannian Manifolds

Description

We provide tools for an isotropic spherical normal (SN) distributions on a (p-1)-sphere in \mathbf{R}^p for sampling, density evaluation, and maximum likelihood estimation of the parameters where the density is defined as

f_SN(x; μ, λ) = \frac{1}{Z(λ)} \exp ≤ft( -\frac{λ}{2} d^2(x,μ) \right)

for location and concentration parameters μ and λ respectively and the normalizing constant Z(λ).

Usage

 1 2 3 4 5 dspnorm(data, mu, lambda, log = FALSE) rspnorm(n, mu, lambda) mle.spnorm(data, method = c("Newton", "Halley", "Optimize", "DE"), ...)

Arguments

 data data vectors in form of either an (n\times p) matrix or a length-n list. See wrap.sphere for descriptions on supported input types. mu a length-p unit-norm vector of location. lambda a concentration parameter that is positive. log a logical; TRUE to return log-density, FALSE for densities without logarithm applied. n the number of samples to be generated. method an algorithm name for concentration parameter estimation.i ... extra parameters for computations, including maxitermaximum number of iterations to be run (default:50). epstolerance level for stopping criterion (default: 1e-5).

Value

dspnorm gives a vector of evaluated densities given samples. rspnorm generates unit-norm vectors in \mathbf{R}^p wrapped in a list. mle.spnorm computes MLEs and returns a list containing estimates of location (mu) and concentration (lambda) parameters.

Examples

 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 # ------------------------------------------------------------------- # Example with Spherical Normal Distribution # # Given a fixed set of parameters, generate samples and acquire MLEs. # Especially, we will see the evolution of estimation accuracy. # ------------------------------------------------------------------- ## DEFAULT PARAMETERS true.mu = c(1,0,0,0,0) true.lbd = 5 ## GENERATE DATA N=1000 big.data = rspnorm(1000, true.mu, true.lbd) ## ITERATE FROM 50 TO 1000 by 10 idseq = seq(from=50, to=1000, by=10) nseq = length(idseq) hist.mu = rep(0, nseq) hist.lbd = rep(0, nseq) for (i in 1:nseq){ small.data = big.data[1:idseq[i]] # data subsetting small.MLE = mle.spnorm(small.data) # compute MLE hist.mu[i] = acos(sum(small.MLE$mu*true.mu)) # difference in mu hist.lbd[i] = small.MLE$lambda } ## VISUALIZE opar <- par(no.readonly=TRUE) par(mfrow=c(1,2)) plot(idseq, hist.mu, "b", pch=19, cex=0.5, main="difference in location") plot(idseq, hist.lbd, "b", pch=19, cex=0.5, main="concentration param") abline(h=true.lbd, lwd=2, col="red") par(opar)

Riemann documentation built on June 20, 2021, 5:07 p.m.