# R/mu_solver.R In SparseMDC: Implementation of SparseMDC Algorithm

```#' Mu Solver
#'
#' Calculates the regularized center values for a cluster. This function
#'  runs inside \code{sparse_mdc}.
#'
#' @param d - The current condition.
#' @param mu - A matrix of regularized mean values.
#' @param v - Relationship matrix, describing relationship between d-1 and d.
#' @param EQ - Equality matrix specifying the number of following terms to
#' which mu_d is equal.
#' @param dim - Total dimensions.
#' @param x - Matrix of mean values.
#' @param nk - Vector with the number of samples in each condition for this
#' cluster.
#' @param p1 - Penalties on mean size.
#' @param p2 - Penalties on mean difference.
#' @param delta Small term to ensure existance of solution.
#'
#' @return A matrix containing the regularized cluster means for each
#' dimension.
#' @export
mu_solver <- function(d, mu, v, EQ,  dim, x, nk, p1, p2, delta){
if ( d == 1){ # Calculate mu_1

mu[,1] <- mu_calc(d = 1, v = v, EQ = EQ, dim = dim, x = x, nk = nk,
p1 = p1, p2 = p2,  delta = delta)[] # Use mu_calc

}
else if (1 < d){
## Calculate mu_{d-1}
v_1 <- v_2 <- v # Create temporary relationship matrices
v_1[,d-1] <- -1 # mu_{d-1} > mu_d
v_2[,d-1] <- 1 # mu_{d-1} < mu_d
EQ_b <- EQ # Create temporary equality matrix
EQ_b[,d-1] <- EQ[,d] + 1 # Add in new equality
# Calculate mu_{d-1}^{1}
fit_1 <- mu_solver(d = d-1, mu = mu, v = v_1, EQ = EQ, dim = dim, x = x,
nk = nk, p1 = p1, p2 = p2, delta = delta)
# Calculate mu_{d-1}^{2}
fit_2 <- mu_solver(d = d-1, mu = mu, v = v_2, EQ = EQ, dim = dim, x = x,
nk = nk, p1 = p1, p2 = p2, delta = delta)
# Calculate mu_{d-1}^{3}
fit_3 <- mu_solver(d = d - 1, mu = mu, v = v, EQ = EQ_b, dim = dim, x = x,
nk = nk, p1 = p1, p2 = p2, delta = delta)

## Calculate mu_d
t_vals <- mu_calc(d = d, v = v, EQ = EQ, dim = dim, x = x, nk = nk,
p1 = p1, p2 = p2, delta = delta)

## Test Relationships between mu_{d-1} & mu_d and assign results
ind_1 <- fit_1[,d-1] > t_vals[] # mu_{d-1}^{1} > mu_d^{1}
mu[ind_1,1:(d-1)] <- fit_1[ind_1,1:(d-1)]
mu[ind_1,d] <- t_vals[][ind_1]

ind_2 <- fit_2[,d-1] < t_vals[] # mu_{d-1}^{2} < mu_d^{2}
mu[ind_2,1:(d-1)] <- fit_2[ind_2,1:(d-1)]
mu[ind_2,d] <- t_vals[][ind_2]

## Calcualte Equal vals
ind_3 <- !ind_1 & !ind_2 # Store indices of equals
if (sum(ind_3) > 0){ # If any terms are equals
mu[ind_3,1:(d-1)] <- fit_3[ind_3,1:(d-1)] # Assign results
mu[ind_3,d] <- mu[ind_3, d-1]
}
}
return(mu)
}
```

## Try the SparseMDC package in your browser

Any scripts or data that you put into this service are public.

SparseMDC documentation built on May 2, 2019, 4:01 a.m.