Skip to content

Commit

Permalink
Version 0.7.0
Browse files Browse the repository at this point in the history
  • Loading branch information
maiermarco committed Nov 25, 2019
1 parent f79e021 commit 07b29b9
Show file tree
Hide file tree
Showing 22 changed files with 241 additions and 132 deletions.
9 changes: 9 additions & 0 deletions .Rbuildignore
Original file line number Diff line number Diff line change
Expand Up @@ -11,3 +11,12 @@ tests\/testthat\/\_res.*\.R$
.*\.o$
.*\.dll$
.*\.bak$
^\.travis\.yml$
^\.codecov\.yml$
^\.covrignore$
.*/\.gitattributes$
.*/\.gitignore$
^README\.Rmd$
^README\.md$
^.*\.Rproj$
^\.Rproj\.user$
1 change: 1 addition & 0 deletions .codecov.yml
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
comment: false
1 change: 1 addition & 0 deletions .covrignore
Original file line number Diff line number Diff line change
@@ -0,0 +1 @@
R/zzz.R
7 changes: 7 additions & 0 deletions .gitignore
Original file line number Diff line number Diff line change
@@ -0,0 +1,7 @@
.Rproj.user
.Rhistory
.RData
.Ruserdata
src/*.o
src/*.so
src/*.dll
32 changes: 27 additions & 5 deletions .travis.yml
Original file line number Diff line number Diff line change
@@ -1,9 +1,31 @@
# R for travis: see documentation at https://docs.travis-ci.com/user/languages/r

language: r
cache: packages

r:
- oldrel
- release
- devel
cache:
packages: true

# addons needed for "rgl":
addons:
apt:
packages:
- freeglut3
- freeglut3-dev
- libglew1.5
- libglew1.5-dev
- libglu1-mesa
- libglu1-mesa-dev
- libgl1-mesa-glx
- libgl1-mesa-dev
- libfreetype6-dev

r_packages:
- covr

matrix:
include:
- r: oldrel
- r: release
after_success:
- Rscript -e 'library(covr); codecov()'
- r: devel
7 changes: 3 additions & 4 deletions DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
Package: DirichletReg
Type: Package
Version: 0.6-4
Date: 2019-10-09
Version: 0.7-0
Date: 2019-11-25
Title: Dirichlet Regression in R
Description: Implements Dirichlet regression models in R.
Authors@R: person(given = c("Marco", "Johannes"),
Expand All @@ -13,8 +13,7 @@ Depends: R (>= 3.0.0), Formula
Imports: stats, graphics, methods, maxLik
Suggests: rgl, knitr, formatR, testthat
License: GPL (>= 2)
URL: http://dirichletreg.r-forge.r-project.org/
http://r-forge.r-project.org/projects/dirichletreg/
URL: https://github.com/maiermarco/DirichletReg
https://CRAN.R-project.org/package=DirichletReg
ByteCompile: yes
LazyLoad: yes
Expand Down
23 changes: 23 additions & 0 deletions DirichletReg.Rproj
Original file line number Diff line number Diff line change
@@ -0,0 +1,23 @@
Version: 1.0

RestoreWorkspace: No
SaveWorkspace: No
AlwaysSaveHistory: No

EnableCodeIndexing: Yes
UseSpacesForTab: Yes
NumSpacesForTab: 2
Encoding: UTF-8

RnwWeave: Sweave
LaTeX: pdfLaTeX

AutoAppendNewline: Yes
StripTrailingWhitespace: Yes

BuildType: Package
PackageUseDevtools: Yes
PackageInstallArgs: --no-multiarch --with-keep.source
PackageBuildArgs: --resave-data --compact-vignettes=gs+qpdf --md5
PackageBuildBinaryArgs: --build
PackageCheckArgs: --as-cran
4 changes: 3 additions & 1 deletion NEWS
Original file line number Diff line number Diff line change
@@ -1,11 +1,13 @@
News for Package 'DirichletReg'

Changes in Version 0.6-4:
Changes in Version 0.7-0:

� Fixed: Cutoff in formula-deparsing removed.

� Removed 'rgl' dependency (package will be loaded if necessary).

� Moved to GitHub.

Changes in Version 0.6-3:

� Instead of producing an error, density functions now return 'NaN'
Expand Down
51 changes: 19 additions & 32 deletions R/DirichReg.R
Original file line number Diff line number Diff line change
Expand Up @@ -69,18 +69,6 @@ if(verbosity > 0){



# if(ifelse(!missing(data), deparse(oformula[[2]]) %in% names(data), FALSE)){
# Y_full <- data[[deparse(oformula[[2]])]] # GET THE FULL DRData OBJECT FOR INFORMATIONS (Y_full)
# } else { # handle responses in .GlobalEnv
# Y_full <- get(deparse(oformula[[2]]), environment(oformula))
# if(missing(data)){
# data[[deparse(oformula[[2]])]] <- Y_full
# } else {
# assign(deparse(oformula[[2]]), Y_full)
# }
# }
# if(class(Y_full) != "DirichletRegData") stop("the response must be prepared by DR_data")

#<<< get Y #####################################################################

repar <- ifelse(model == "common", FALSE, TRUE)
Expand Down Expand Up @@ -114,8 +102,8 @@ if(verbosity > 0){
y_in <- (seq_len(ncol(Y)))[sub.comp]
y_out <- (seq_len(ncol(Y)))[-sub.comp]
y_in_labels <- colnames(Y)[y_in]
y_out_labels <- paste(colnames(Y)[y_out], sep="", collapse=" + ")
Y <- cbind(rowSums(Y[,y_out]), Y[,y_in])
y_out_labels <- paste(colnames(Y)[y_out], sep = "", collapse = " + ")
Y <- cbind(rowSums(Y[, y_out]), Y[, y_in])
colnames(Y) <- c(y_out_labels, y_in_labels)
}

Expand All @@ -137,12 +125,12 @@ if(verbosity > 0){


if(!repar){
X.mats <- lapply(seq_len(ncol(Y)), function(i){ model.matrix(terms(formula, data=data, rhs=i), mf) })
X.mats <- lapply(seq_len(ncol(Y)), function(i){ model.matrix(terms(formula, data = data, rhs = i), mf) })
Z.mat <- NULL
n.vars <- unlist(lapply(X.mats, ncol))
} else {
X.mats <- lapply(seq_len(ncol(Y)), function(i) model.matrix(terms(formula, data=data, rhs=1), mf) )
Z.mat <- model.matrix(terms(formula, data=data, rhs=2), mf)
X.mats <- lapply(seq_len(ncol(Y)), function(i) model.matrix(terms(formula, data = data, rhs = 1), mf) )
Z.mat <- model.matrix(terms(formula, data = data, rhs = 2), mf)
n.vars <- c(unlist(lapply(X.mats, ncol))[-1], ncol(Z.mat))
}

Expand Down Expand Up @@ -179,9 +167,9 @@ if(verbosity > 0){

# compute starting values
if(is.null(control$sv)){
starting.vals <- get_starting_values(Y=Y_fit, X.mats=X_fit,
Z.mat={if(repar) as.matrix(Z.mat) else Z.mat},
repar=repar, base=base, weights=weights) * if(repar){ 1 } else { 1/n.dim }
starting.vals <- get_starting_values(Y = Y_fit, X.mats = X_fit,
Z.mat = {if(repar) as.matrix(Z.mat) else Z.mat},
repar = repar, base = base, weights = weights) * if(repar){ 1 } else { 1/n.dim }
} else {
if(length(control$sv) != n.vars) stop("wrong number of starting values supplied.")
starting.vals <- control$sv
Expand Down Expand Up @@ -213,7 +201,7 @@ if(verbosity > 0){
coefs <- fit.res$estimate

if(repar){
names(coefs) <- unlist(as.vector(c(rep(colnames(X.mats[[1]]),n.dim-1),colnames(Z.mat))))
names(coefs) <- unlist(as.vector(c(rep(colnames(X.mats[[1]]), n.dim-1), colnames(Z.mat))))
} else {
names(coefs) <- unlist(lapply(X.mats, colnames))
}
Expand All @@ -223,9 +211,9 @@ if(verbosity > 0){
if(repar){

B <- matrix(0.0, nrow = n.vars[1L], ncol = n.dim)
B[cbind(rep(seq_len(n.vars[1L]), (n.dim-1L)), rep(seq_len(n.dim)[-base], each=n.vars[1]))] <- coefs[1:((n.dim-1)*n.vars[1])]
B[cbind(rep(seq_len(n.vars[1L]), (n.dim-1L)), rep(seq_len(n.dim)[-base], each = n.vars[1]))] <- coefs[1:((n.dim-1)*n.vars[1])]

g <- matrix(coefs[((n.dim-1)*n.vars[1]+1):length(coefs)], ncol=1)
g <- matrix(coefs[((n.dim-1)*n.vars[1]+1):length(coefs)], ncol = 1)

XB <- exp(apply(B, 2L, function(b){ as.matrix(X.mats[[1L]]) %*% b }))
MU <- apply(XB, 2L, function(x){ x /rowSums(XB) })
Expand All @@ -236,9 +224,9 @@ if(verbosity > 0){

} else {

B <- sapply(seq_len(n.dim), function(i){ coefs[(cumsum(c(0,n.vars))[i]+1) : cumsum(n.vars)[i]] }, simplify = FALSE)
B <- sapply(seq_len(n.dim), function(i){ coefs[(cumsum(c(0, n.vars))[i]+1) : cumsum(n.vars)[i]] }, simplify = FALSE)

ALPHA <- sapply(seq_len(n.dim), function(i){ exp(as.matrix(X.mats[[i]]) %*% matrix(B[[i]], ncol=1)) })
ALPHA <- sapply(seq_len(n.dim), function(i){ exp(as.matrix(X.mats[[i]]) %*% matrix(B[[i]], ncol = 1)) })

PHI <- rowSums(ALPHA)
MU <- apply(ALPHA, 2L, "/", PHI)
Expand All @@ -251,21 +239,21 @@ if(verbosity > 0){
hessian <- fit.res$hessian

vcov <- tryCatch(solve(-fit.res$hessian),
error=function(x){ return(matrix(NA, nrow=nrow(hessian), ncol=ncol(hessian))) },
silent=TRUE)
error = function(x){ return(matrix(NA, nrow = nrow(hessian), ncol = ncol(hessian))) },
silent = TRUE)

if(!repar){ ## COMMON
coefnames <- apply(cbind(rep(varnames, n.vars), unlist(lapply(X.mats, colnames))), 1, paste, collapse=":")
coefnames <- apply(cbind(rep(varnames, n.vars), unlist(lapply(X.mats, colnames))), 1, paste, collapse = ":")
} else { ## ALTERNATIVE
coefnames <- apply(cbind(rep(c(varnames[-base], "(phi)"), n.vars), c(unlist(lapply(X.mats, colnames)[-base]), colnames(Z.mat))), 1, paste, collapse=":")
coefnames <- apply(cbind(rep(c(varnames[-base], "(phi)"), n.vars), c(unlist(lapply(X.mats, colnames)[-base]), colnames(Z.mat))), 1, paste, collapse = ":")
}

dimnames(hessian) <- list(coefnames, coefnames)
dimnames(vcov) <- list(coefnames, coefnames)
shortnames <- names(coefs)
names(coefs) <- coefnames

se <- if(!any(is.na(vcov))) sqrt(diag(vcov)) else rep(NA,length(coefs))
se <- if(!any(is.na(vcov))) sqrt(diag(vcov)) else rep(NA, length(coefs))

res <- structure(list(
call = this.call,
Expand All @@ -287,7 +275,7 @@ if(verbosity > 0){
npar = length(coefs),
coefficients = coefs,
coefnames = shortnames,
fitted.values = list(mu=MU,phi=PHI,alpha=ALPHA),
fitted.values = list(mu = MU, phi = PHI, alpha = ALPHA),
logLik = fit.res$maximum,
vcov = vcov,
hessian = hessian,
Expand All @@ -311,5 +299,4 @@ if(verbosity > 0){

return(res)


}
41 changes: 15 additions & 26 deletions R/Dirichlet.R
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
rdirichlet <- function(n, # single integer specifying the sample size
alpha # can be a single vector or a matrix (with exactly n rows)
){

rdirichlet <- function(
n # single integer specifying the sample size
, alpha # can be a single vector or a matrix (with exactly n rows)
){
# check if the sample size is an integer > 0
if( ((n %% 1) != 0) | (n <= 0)) stop("n must be an integer > 0")
# check if any value in alpha is <= 0
Expand All @@ -10,32 +10,20 @@ rdirichlet <- function(n, # single integer specifying the sample size
.vec <- is.vector(alpha)
.mat <- is.matrix(alpha)

if(!.vec & !.mat){ # alpha is neither a vector nor a matrix

if(!.vec && !.mat){ # alpha is neither a vector nor a matrix
stop("alpha must be a vector or a matrix")

} else if(.vec & !.mat){ # alpha is a vector

# dims <- length(alpha)
# G <- matrix(rgamma(n*dims, alpha, 1), byrow=TRUE, ncol=dims)
# X <- G / rowSums(G)
} else if(.vec && !.mat){ # alpha is a vector
X <- .Call("rdirichlet_vector", n, alpha)

} else { # alpha is a matrix

# dims <- ncol(alpha)
if(n != nrow(alpha)) stop("when alpha is a matrix, the number of its rows must be equal to n")
#
# G <- matrix(rgamma(n*dims, as.vector(alpha), 1), nrow=n)
# X <- G / rowSums(G)
} else { # alpha is a matrix
if(n != nrow(alpha)){
stop("when alpha is a matrix, the number of its rows must be equal to n")
}
X <- .Call("rdirichlet_matrix", n, alpha, dim(alpha))

}

return(X)


}###END OF rdirichlet
}



Expand Down Expand Up @@ -65,20 +53,21 @@ ddirichlet <- function(x, alpha, log = FALSE, sum.up = FALSE){



ddirichlet_R <- function(x, alpha, log = FALSE, sum.up = FALSE){ # PURE R VERSION
# PURE R VERSION
ddirichlet_R <- function(x, alpha, log = FALSE, sum.up = FALSE){
# some checking!
if(is.null(dim(x))) stop("x must be a matrix")
if(is.vector(alpha)){
if(ncol(x) != length(alpha)) stop("alpha must be a vector/matrix fitting to the data in x")
alpha <- matrix(rep(alpha,nrow(x)),nrow(x),byrow=T)
alpha <- matrix(rep(alpha, nrow(x)), nrow(x), byrow = TRUE)
}
if(any(dim(alpha) != dim(x))) stop("check if x and alpha are correctly specified")
if(any(alpha <= 0)){
warning("all values in alpha must be > 0")
if(sum.up) return(NaN) else return(rep(NaN, nrow(x)))
}

res <- lgamma(rowSums(alpha)) - rowSums(lgamma(alpha)) + rowSums((alpha-1)*log(x))
res <- lgamma(rowSums(alpha)) - rowSums(lgamma(alpha)) + rowSums((alpha-1.0)*log(x))

if(sum.up){
if(log) return(sum(res)) else return(exp(sum(res)))
Expand Down
Loading

0 comments on commit 07b29b9

Please sign in to comment.