Nothing
#=================================================================================#
# UNIFORM RANDOM MATRICES
#=================================================================================#
#' @title Generate a uniform random matrix
#' @description Uniform random matrices are matrices with uniformly distributed entries. They are an elementary type of random matrix.
#'
#' @param N number of dimensions of the square matrix
#' @param min minimum of the uniform distribution to be sampled from
#' @param max maximum of the uniform distribution to be sampled from
#' @param symm indicates whether the matrix should be symmetric (equal to its transpose).
#' @param cplx indicates whether the matrix should have complex entries.
#' @param herm indicates whether the matrix should be hermitian (equal to its conjugate transpose).
#' Reserved for when cplx = TRUE, otherwise use symm = TRUE.
#'
#' @return A random matrix with uniformly distributed entries.
#'
#' @examples
#' # Unif(1,2) distributed matrix
#' P <- RM_unif(N = 3, min = 1, max = 2)
#'
#' # Unif(0,5) distributed matrix with real symmetric entries
#' P <- RM_unif(N = 7, min = 0, max = 5, symm = TRUE)
#'
#' # Unif(0,1) distributed matrix with complex entries
#' Q <- RM_unif(N = 7, min = 0, max = 1, cplx = TRUE)
#'
#' # Unif(2,10) distributed matrix with hermitian complex entries
#' Q <- RM_unif(N = 5, min = 2, max = 10, cplx = TRUE, herm = TRUE)
#'
RM_unif <- function(N, min, max, symm = FALSE, cplx = FALSE, herm = FALSE){
# Create [N x N] matrix with uniformly distributed entries
P <- matrix(runif(N^2, min, max), nrow = N)
# Make symmetric/hermitian if prompted
if(symm || herm){P <- .makeHermitian(P)}
# Returns a matrix with complex (and hermitian) entries if prompted
if(cplx){
Im_P <- (1i * RM_unif(N, min, max)) # Recursively add imaginary components as 1i * instance of real-valued matrix.
if(herm){P <- P + .makeHermitian(Im_P)} # Make imaginary part hermitian if prompted
else{P <- P + Im_P}
}
P # Return the matrix
}
#=================================================================================#
# NORMAL RANDOM MATRICES
#=================================================================================#
#' @title Generate a normal random matrix
#' @description Normal random matrices are matrices with normally distributed entries. These matrices
#' are extensively studied in random matrix theory.
#'
#' @param N number of dimensions of the square matrix
#' @param mean mean of the normal distribution of entries
#' @param sd standard deviation of the normal distribution of entries
#' @param symm indicates whether the matrix should be symmetric (equal to its transpose).
#' Reserved for when cplx = FALSE, otherwise use herm = TRUE.
#' @param cplx indicates whether the matrix should have complex entries.
#' @param herm indicates whether the matrix should be hermitian (equal to its conjugate transpose).
#' Reserved for when cplx = TRUE, otherwise use symm = TRUE.
#'
#' @return A random matrix with normally distributed entries.
#' @examples
#' # N(1,2) distributed matrix
#' P <- RM_norm(N = 3, mean = 1, sd = 2)
#'
#' # N(0,5) distributed matrix with real symmetric entries
#' P <- RM_norm(N = 7, sd = 5, symm = TRUE)
#'
#' # 7x7 standard normal matrix with complex entries
#' Q <- RM_norm(N = 7, cplx = TRUE)
#'
#' # N(2,1) distributed matrix with hermitian complex entries
#' Q <- RM_norm(N = 5, mean = 2, cplx = TRUE, herm = TRUE)
#'
RM_norm <- function(N, mean = 0, sd = 1, symm = FALSE, cplx = FALSE, herm = FALSE){
# Create [N x N] matrix with normally distributed entries
P <- matrix(rnorm(N^2, mean, sd), nrow = N)
# Make symmetric/hermitian if prompted
if(symm || herm){P <- .makeHermitian(P)}
# Returns a matrix with complex (and hermitian) entries if prompted
if(cplx){
Im_P <- (1i * RM_norm(N, mean, sd)) # Recursively add imaginary components as 1i * instance of real-valued matrix.
if(herm){P <- P + .makeHermitian(Im_P)} # Make imaginary part hermitian if prompted
else{P <- P + Im_P}
}
P # Return the matrix
}
#=================================================================================#
#' @title Generate a Hermite \eqn{\beta}-matrix
#' @description Hermite-\eqn{\beta} ensemble matrices are matrices with normal entries and beta real number components.
#' Using Ioana Dumitriu's tridiagonal matrix model, this function is an implementation of the generalized, but not necessarily invariant,
#' beta ensembles for \eqn{\beta} > 0.
#'
#' @param N number of dimensions of the square matrix
#' @param beta the value of the beta parameter for the beta ensemble
#'
#' @return A random Hermite beta matrix with any integer parameter beta
#'
#' @examples
#' # Generate a 3x3 random beta matrix with beta = 4
#' P <- RM_beta(N = 3, beta = 4)
#'
#' # Generate a 10x10 random beta matrix with beta = 25
#' P <- RM_beta(N = 10, beta = 25)
#'
RM_beta <- function(N, beta){
# Set the diagonal ~ N(0,2)
P <- diag(rnorm(n = N, mean = 0, sd = sqrt(2)))
# Get degrees of freedom sequence for offdigonal
df_seq <- beta * (N - seq(1, N-1))
# Set the off-1 diagonals as chi squared variables with df(beta_i)
P[row(P) - col(P) == 1] <- P[row(P) - col(P) == -1] <- sqrt(rchisq(N-1, df_seq))
# Rescale the entries by 1/sqrt(2)
P <- P/sqrt(2)
# Return the beta matrix
P
}
#=================================================================================#
#' @title Generate a tridiagonal matrix with normal entries
#'
#' @param N number of dimensions of the square matrix
#' @param symm indicates whether the matrix should be symmetric; equal to its transpose.
#'
#' @return A random tridiagonal matrix with N(0,2) diagonal and N(0,1) band.
#'
#' @examples
#' # Generate a 3x3 standard normal tridiagonal matrix
#' P <- RM_trid(N = 3)
#'
#' # Symmetric tridiagonal matrix
#' P <- RM_trid(N = 9, symm = TRUE)
#'
RM_trid <- function(N, symm = FALSE){
diagonal <- rnorm(n = N, 0, 2)
P <- diag(diagonal)
P[row(P) - col(P) == 1] <- P[row(P) - col(P) == -1] <- rnorm(n = N, 0, 1)
P # Return the matrix
}
#=================================================================================#
# STOCHASTIC RANDOM MATRICES
#=================================================================================#
#' @title Generate a random stochastic matrix
#' @description A (row-)stochastic matrix is a matrix whose rums sum to 1. There is a natural one-to-one corrospondence between
#' stochastic matrices and Markov Chains; this is so when its i,j entry represent the transition probability from state i to state j.
#'
#' @param N number of dimensions of the square matrix
#' @param symm indicates whether the matrix should be symmetric; equal to its transpose.
#' @param sparsity indicates whether the matrix should add some arbitrary sparsity (zeros)
#'
#' @return A random stochastic matrix.
#'
#' @examples
#' P <- RM_stoch(N = 3)
#' P <- RM_stoch(N = 9, sparsity = TRUE)
#' Q <- RM_stoch(N = 9, symm = TRUE)
#' Q <- RM_stoch(N = 9, symm = TRUE, sparsity = TRUE)
#'
RM_stoch <- function(N, symm = FALSE, sparsity = FALSE){
# Choose row function depending on sparsity argument
if(sparsity){row_fxn <- .stoch_row_zeros} else {row_fxn <- .stoch_row}
# Generate the [N x N] stochastic matrix stacking N stochastic rows
P <- do.call("rbind", lapply(X = rep(N, N), FUN = row_fxn))
# Make symmetric (if prompted)
if(symm){ P <- .makeStochSymm(P) }
# Return the matrix
P
}
#=================================================================================#
#' @title Generate a random stochastic matrix for a walk on an Erdos-Renyi graph
#' @description An Erdos-Renyi Graph is a graph whose edges are connected ~ Bern(p).
#' Hence, its transition matrix will have nonzero entries with that probability.
#' So, we can alternatively think of the transition matrix for such walk as a stochastic matrix with parameterized sparsity.
#'
#' @param N number of dimensions of the square matrix
#' @param p the probability two vertices are connected in an Erdos-Renyi graph.
#'
#' @return A random stochastic matrix corrosponding to a walk on an Erdos-Renyi graph with probability p.
#'
#' @examples
#' # Very sparse graph
#' P <- RM_erdos(N = 3, p = 0.2)
#'
#' # Slightly sparse graph
#' P <- RM_erdos(N = 9, p = 0.6)
#'
#' # Completely connected graph
#' P <- RM_erdos(N = 5, p = 1)
#'
RM_erdos <- function(N, p){
# Generate an [N x N] Erdos-Renyi stochastic matrix by stacking N p-stochastic rows
P <- do.call("rbind", lapply(X = rep(N, N), FUN = .stoch_row_erdos, p = p))
# Return the Erdos-Renyi transition matrix
P
}
#=================================================================================#
# STOCHASTIC ROW FUNCTIONS
#=================================================================================#
# Generates stochastic rows of length N
.stoch_row <- function(N){
# Sample a vector of probabilities
row <- runif(n = N, min = 0, max = 1)
# Return the normalized row (sums to one)
row / sum(row)
}
#=================================================================================#
# Generates same rows as in .stoch_row(N), but with randomly introduced sparsity
.stoch_row_zeros <- function(N){
# Sample a vector of probabilities
row <- runif(n = N, min = 0, max = 1)
# Sample a vertex degree of at least one (as to ensure row is stochastic)
degree_vertex <- sample(x = 1:(N-1), size = 1)
# Sever a random selection of edges to set the vertex degree
row[sample(1:N, size = N - degree_vertex)] <- 0
# Return normalized row
row / sum(row)
}
#=================================================================================#
# Generates a stochastic row with parameterized sparsity of p
.stoch_row_erdos <- function(N, p){
# Sample a vector of probabilities
row <- runif(n = N, min = 0, max = 1)
# Sample the vertex degree so that it is ~ Bin(n,p)
degree_vertex <- rbinom(n = 1, size = N, prob = 1 - p)
# Sever a random selection of edges to set the vertex degree
row[sample(1:N, degree_vertex)] <- 0
# Return normalized row only if non-zero (cannot divide by 0)
if(sum(row) != 0){
row / sum(row)
} else{
.stoch_row_erdos(N, p) # Otherwise, try again
}
}
#=========================================================================#
# HELPER FUNCTIONS
#=========================================================================#
# Returns a Hermitian version of a matrix by manual assignment
.makeHermitian <- function(P){
for(i in 1:nrow(P)){
for(j in 1:ncol(P)){
# Select the entries in the upper triangle (i < j)
if(i < j){
# Make the upper triangle equal to the conjugate transpose of the lower triangle
P[i,j] <- Conj(P[j,i])
}
}
}
# Return the Hermitian Matrix
P
}
#=================================================================================#
# Take a stochastic matrix and make it symmetric
.makeStochSymm <- function(P_){
# Get parameters
N <- nrow(P_)
# Make lower and upper triangles equal to each other's conjugate transpose
P <- .makeHermitian(P_)
# Nullify diagonal
diag(P) <- rep(0, N)
# Normalize rows
for(i in 1:N){P[i, ] <- P[i, ]/sum(P[i, ])}
# Set diagonal to diffrence between 1 and sum of offdigonal entries to make rows sum to 1
diag <- vector(mode = "numeric", length = N)
for(i in 1:N){diag[i] <- (1 - sum(.offdiagonalEntries(row = P[i, ], row_index = i)))}
diag(P) <- diag
# Return the symmetric stochastic matrix
return(P)
}
#=================================================================================#
# Return the off-diagonal entries of row i
.offdiagonalEntries <- function(row, row_index){
row[which(1:length(row) != row_index)]
}
#=================================================================================#
# RANDOM MATRIX ENSEMBLES
#=================================================================================#
#' @title Generate an ensemble of normal random matrices
#' @description Given the same arguments as RM_norm, this function returns an ensemble of random normal matrices.
#' While random matrices usually do not exude unique properties on their own, they do indeed have
#' deterministic properties at the ensemble level in terms of their spectral statistics.
#'
#' @inheritParams RM_unif
#' @param ... any default-valued parameters taken as arguments by RM_norm()
#' @param size the size of the ensemble (i.e. number of matrices)
#'
#' @return An ensemble (list) of normal matrices as specified by the matrix arguments.
#'
#' @examples
#' # Generate an ensemble of standard normal 3x3 matrices of size 20
#' ensemble <- RME_norm(N = 3, size = 20)
#'
RME_unif <- function(N, min, max, ..., size){lapply(X = rep(N, size), FUN = RM_unif, min, max, ...)}
#=================================================================================#
#' @title Generate an ensemble of normal random matrices
#' @description Given the same arguments as RM_norm, this function returns an ensemble of random normal matrices.
#' While random matrices usually do not exude unique properties on their own, they do indeed have
#' deterministic properties at the ensemble level in terms of their spectral statistics.
#'
#' @inheritParams RM_norm
#' @param ... any default-valued parameters taken as arguments by RM_norm()
#' @param size the size of the ensemble (i.e. number of matrices)
#'
#' @return An ensemble (list) of normal matrices as specified by the matrix arguments.
#'
#' @examples
#' # Generate an ensemble of standard normal 3x3 matrices of size 20
#' ensemble <- RME_norm(N = 3, size = 20)
#'
RME_norm <- function(N, mean = 0, sd = 1, ..., size){lapply(X = rep(N, size), FUN = RM_norm, mean, sd, ...)}
#=================================================================================#
#' @title Generate an ensemble of random beta matrices
#' @description Given the same arguments as RM_norm, this function returns an ensemble of that particular class of matrix.
#' While random matrices usually do not exude unique properties on their own, they do indeed have
#' deterministic properties at the ensemble level in terms of their spectral statistics.
#'
#' @inheritParams RM_beta
#' @param size the size of the ensemble (i.e. number of matrices)
#'
#' @return An ensemble (list) of beta matrices as specified by the matrix arguments.
#'
#' @examples
#' # Generate an ensemble of 10x10 beta matrices with beta = 4 of size 100.
#' ensemble <- RME_beta(N = 10, beta = 4, size = 100)
#'
RME_beta <- function(N, beta, size){lapply(X = rep(N, size), FUN = RM_beta, beta)}
#=================================================================================#
#' @title Generate an ensemble of stochastic matrices
#' @description Given the same arguments as RM_stoch, this function returns an ensemble of random stochastic matrices.
#' While random matrices usually do not exude unique properties on their own, they do indeed have
#' deterministic properties at the ensemble level in terms of their spectral statistics.
#'
#' @inheritParams RM_stoch
#' @param ... pass any default-valued parameters taken as arguments by RM_stoch()
#' @param size the size of the ensemble (i.e. number of matrices)
#'
#' @return An ensemble (list) of stochastic matrices as specified by the matrix arguments.
#'
#' @examples
#' # Generate an ensemble of random 5x5 transition matrices of size 20.
#' ensemble <- RME_stoch(N = 5, size = 20)
#'
#' # Generate an ensemble of symmetric random 5x5 transition matrices of size 20.
#' ensemble <- RME_stoch(N = 5, symm = TRUE, size = 20)
#'
RME_stoch <- function(N, ..., size){lapply(X = rep(N, size), FUN = RM_stoch, ...)}
#=================================================================================#
#' @title Generate an ensemble of Erdos-Renyi transition matrices
#' @description Given the same arguments as RM_norm, this function returns an ensemble of random Erdos-Renyi stochastic matrices.
#' While random matrices usually do not exude unique properties on their own, they do indeed have
#' deterministic properties at the ensemble level in terms of their spectral statistics.
#'
#' @inheritParams RM_erdos
#' @param size the size of the ensemble (i.e. number of matrices)
#'
#' @return An ensemble (list) of Erdos-Renyi transition matrices as specified by the matrix arguments.
#'
#' @examples
#' # Generate an ensemble of 10x10 Erdos-Renyi transition matrices of size 50 with p = 0.7
#' ensemble <- RME_erdos(N = 10, p = 0.7, size = 50)
#'
RME_erdos <- function(N, p, size){lapply(X = rep(N, size), FUN = RM_erdos, p)}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.