Skip to content

Commit

Permalink
fix dependency on VGAM gaussianff, which was removed in 1.06
Browse files Browse the repository at this point in the history
  • Loading branch information
ctrapnell committed Oct 17, 2018
1 parent 2c9def9 commit 3b09319
Show file tree
Hide file tree
Showing 11 changed files with 12 additions and 161 deletions.
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ Depends:
Matrix (>= 1.2-6),
Biobase,
ggplot2 (>= 1.0.0),
VGAM (>= 1.0-1),
VGAM (>= 1.0-6),
DDRTree (>= 0.1.4),
Imports:
parallel,
Expand Down
2 changes: 0 additions & 2 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -58,7 +58,6 @@ export(reducedDimS)
export(reducedDimW)
export(relative2abs)
export(responseMatrix)
export(selectNegentropyGenes)
export(selectTopMarkers)
export(setOrderingFilter)
export(vstExprs)
Expand Down Expand Up @@ -172,7 +171,6 @@ importFrom(stats,optim)
importFrom(stats,p.adjust)
importFrom(stats,prcomp)
importFrom(stats,qnorm)
importFrom(stats,quantile)
importFrom(stats,rnorm)
importFrom(stats,sd)
importFrom(stats,setNames)
Expand Down
2 changes: 1 addition & 1 deletion R/RcppExports.R
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,6 @@
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393

jaccard_coeff <- function(R_idx, R_weight) {
.Call('monocle_jaccard_coeff', PACKAGE = 'monocle', R_idx, R_weight)
.Call('_monocle_jaccard_coeff', PACKAGE = 'monocle', R_idx, R_weight)
}

4 changes: 2 additions & 2 deletions R/cds_conversion.R
Original file line number Diff line number Diff line change
Expand Up @@ -149,7 +149,7 @@ importCDS <- function(otherCDS, import_all = FALSE) {
if(all(data == floor(data))) {
expressionFamily <- negbinomial.size()
} else if(any(data < 0)){
expressionFamily <- gaussianff()
expressionFamily <- uninormal()
} else {
expressionFamily <- tobit()
}
Expand Down Expand Up @@ -232,7 +232,7 @@ importCDS <- function(otherCDS, import_all = FALSE) {
if(all(data == floor(data))) {
expressionFamily <- negbinomial.size()
} else if(any(data < 0)){
expressionFamily <- gaussianff()
expressionFamily <- uninormal()
} else {
expressionFamily <- tobit()
}
Expand Down
2 changes: 1 addition & 1 deletion R/differential_expression.R
Original file line number Diff line number Diff line change
Expand Up @@ -42,7 +42,7 @@ diff_test_helper <- function(x,
expressionFamily <- negbinomial.size(size=1/disp_guess)
}
}
}else if (expressionFamily@vfamily %in% c("gaussianff", "uninormal")){
}else if (expressionFamily@vfamily %in% c("uninormal")){
f_expression <- x
}else if (expressionFamily@vfamily %in% c("binomialff")){
f_expression <- x
Expand Down
27 changes: 3 additions & 24 deletions R/expr_models.R
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ fit_model_helper <- function(x,
}
}
}
else if (expressionFamily@vfamily %in% c("gaussianff", "uninormal", "binomialff")) {
else if (expressionFamily@vfamily %in% c("uninormal", "binomialff")) {
f_expression <- x
}
else {
Expand All @@ -64,7 +64,7 @@ fit_model_helper <- function(x,
if (expressionFamily@vfamily %in% c("negbinomial", "negbinomial.size")){
disp_guess <- calculate_NB_dispersion_hint(disp_func, round(orig_x), expr_selection_func = max)
backup_expression_family <- negbinomial()
}else if (expressionFamily@vfamily %in% c("gaussianff", "uninormal")){
}else if (expressionFamily@vfamily %in% c("uninormal")){
backup_expression_family <- NULL
}else if (expressionFamily@vfamily %in% c("binomialff")){
backup_expression_family <- NULL
Expand Down Expand Up @@ -152,7 +152,7 @@ responseMatrix <- function(models, newdata = NULL, response_type="response", cor
if (is.null(x)) { NA } else {
if (x@family@vfamily %in% c("negbinomial", "negbinomial.size")) {
predict(x, newdata = newdata, type = response_type)
} else if (x@family@vfamily %in% c("gaussianff")) {
} else if (x@family@vfamily %in% c("uninormal")) {
predict(x, newdata = newdata, type = response_type)
}
else {
Expand Down Expand Up @@ -403,27 +403,6 @@ parametricDispersionFit <- function( disp_table, verbose = FALSE, initial_coefs=
list(fit, coefs)
}

# parametricDispersionFit <- function( means, disps )
# {
# coefs <- c( 1e-6, 1 )
# iter <- 0
#
# residuals <- disps / ( coefs[1] + coefs[2] / means )
# good <- which( (residuals > 1e-4) & (residuals < 10000) )
# fit <- vglm( log(disps[good]) ~ log(means[good]), family=gaussianff())
# oldcoefs <- coefs
# coefs <- coefficients(fit)
#
# iter <- iter + 1
# print(coefs)
# names( coefs ) <- c( "asymptDisp", "extraPois" )
# ans <- function( q )
# exp(coefs[1] + coefs[2] * log(q))
# #ans
# coefs
# }


#' Return a variance-stabilized matrix of expression values
#'
#' @description This function was taken from the DESeq package (Anders and Huber) and modified
Expand Down
2 changes: 1 addition & 1 deletion R/order_cells.R
Original file line number Diff line number Diff line change
Expand Up @@ -1263,7 +1263,7 @@ normalize_expr_data <- function(cds,
}else{
stop("Error: the only normalization methods supported with Tobit-distributed (e.g. FPKM/TPM) data are 'log' (recommended) or 'none'")
}
}else if (cds@expressionFamily@vfamily == "gaussianff") {
}else if (cds@expressionFamily@vfamily == "uninormal") {
if (norm_method == "none"){
FM <- FM + pseudo_expr
}else{
Expand Down
90 changes: 0 additions & 90 deletions R/utils.R
Original file line number Diff line number Diff line change
Expand Up @@ -206,96 +206,6 @@ smartEsApply <- function(X, MARGIN, FUN, convert_to_dense, ...) {
}


####
#' Filter genes with extremely high or low negentropy
#'
#' @description Examines all the genes in the CellDataSet passed in and removes
#' all the genes that contain extremely high or low negentropies. You can specify
#' which genes to filter out based on the boundaries you can set for expression levels
#' and the boundaries you set for which centile to include. the function "dispersionTable"
#' is a better form of this function.
#'
#' @param cds a CellDataSet object upon which to perform this operation
#' @param lower_negentropy_bound the centile below which to exclude to genes
#' @param upper_negentropy_bound the centile above which to exclude to genes
#' @param expression_lower_thresh the expression level below which to exclude genes used to determine negentropy
#' @param expression_upper_thresh the expression level above which to exclude genes used to determine negentropy
#' @return a vector of gene names
#' @importFrom stats quantile
#' @export
#' @examples
#' \dontrun{
#' reasonableNegentropy <- selectNegentropyGenes(HSMM, "07%", "95%", 1, 100)
#' }
selectNegentropyGenes <- function(cds, lower_negentropy_bound="0%",
upper_negentropy_bound="99%",
expression_lower_thresh=0.1,
expression_upper_thresh=Inf){
.Deprecated("dispersionTable")
log_expression <- NULL

FM <- exprs(cds)
if (cds@expressionFamily@vfamily %in% c("negbinomial", "negbinomial.size"))
{
expression_lower_thresh <- expression_lower_thresh / colSums(FM)
expression_upper_thresh <- expression_upper_thresh / colSums(FM)
FM <- Matrix::t(Matrix::t(FM)/colSums(FM))
}

negentropy_exp <- apply(FM,1,function(x) {

expression <- x[x > expression_lower_thresh]
expression <- log2(x);
expression <- expression[is.finite(expression)]

if (length(expression)){
expression <- scale(expression)
mean(-exp(-(expression^2)/2))^2
}else{
0
}
}

)


means <- apply(FM,1,function(x) {
expression <- x[x > expression_lower_thresh]
expression <- log2(x);
expression[is.finite(expression) == FALSE] <- NA;

if (length(expression)){
mean(expression, na.rm=T)
}else{
NA
}
}
)
ordering_df <- data.frame(log_expression = means, negentropy = negentropy_exp)

negentropy <- NULL
log_express <- NULL
negentropy_residual <- NULL

ordering_df <- subset(ordering_df,
is.na(log_expression) == FALSE &
is.nan(log_expression) == FALSE &
is.na(negentropy) == FALSE &
is.nan(negentropy) == FALSE)
negentropy_fit <- vglm(negentropy~sm.ns(log_expression, df=4),data=ordering_df, family=VGAM::gaussianff())
ordering_df$negentropy_response <- predict(negentropy_fit, newdata=ordering_df, type="response")
ordering_df$negentropy_residual <- ordering_df$negentropy - ordering_df$negentropy_response
lower_negentropy_thresh <- quantile(ordering_df$negentropy_residual, probs=seq(0,1,0.01), na.rm=T)[lower_negentropy_bound]
upper_negentropy_thresh <- quantile(ordering_df$negentropy_residual, probs=seq(0,1,0.01), na.rm=T)[upper_negentropy_bound]

ordering_genes <- row.names(subset(ordering_df,
negentropy_residual >= lower_negentropy_thresh &
negentropy_residual <= upper_negentropy_thresh))
ordering_genes
}



#' Retrieve a table of values specifying the mean-variance relationship
#'
#' Calling estimateDispersions computes a smooth function describing how variance
Expand Down
2 changes: 1 addition & 1 deletion actual_vignette_holder/monocle-vignette-original.Rnw
Original file line number Diff line number Diff line change
Expand Up @@ -422,7 +422,7 @@ Monocle works well with both relative expression data and count-based measures (
\hline
\noalign{\smallskip}

\Rfunction{gaussianff()} & log-transformed FPKM/TPMs, Ct values from single-cell qPCR & If you want to use Monocle on data you have already transformed to be normally distributed, you can use this function, though some Monocle features may not work well.
\Rfunction{uninormal()} & log-transformed FPKM/TPMs, Ct values from single-cell qPCR & If you want to use Monocle on data you have already transformed to be normally distributed, you can use this function, though some Monocle features may not work well.
\end{tabular}

\textbf{Using the wrong expressionFamily for your data will lead to bad results}, errors from Monocle, or both. However, if you have FPKM/TPM data, you can still use negative binomial if you first convert your relative expression values to transcript counts using \Rfunction{relative2abs}. This often leads to much more accurate results than using \Rfunction{tobit()}. See \ref{conv_rel2abs} for details.
Expand Down
36 changes: 0 additions & 36 deletions man/selectNegentropyGenes.Rd

This file was deleted.

4 changes: 2 additions & 2 deletions src/RcppExports.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -7,7 +7,7 @@ using namespace Rcpp;

// jaccard_coeff
NumericMatrix jaccard_coeff(SEXP R_idx, SEXP R_weight);
RcppExport SEXP monocle_jaccard_coeff(SEXP R_idxSEXP, SEXP R_weightSEXP) {
RcppExport SEXP _monocle_jaccard_coeff(SEXP R_idxSEXP, SEXP R_weightSEXP) {
BEGIN_RCPP
Rcpp::RObject rcpp_result_gen;
Rcpp::RNGScope rcpp_rngScope_gen;
Expand All @@ -19,7 +19,7 @@ END_RCPP
}

static const R_CallMethodDef CallEntries[] = {
{"monocle_jaccard_coeff", (DL_FUNC) &monocle_jaccard_coeff, 2},
{"_monocle_jaccard_coeff", (DL_FUNC) &_monocle_jaccard_coeff, 2},
{NULL, NULL, 0}
};

Expand Down

0 comments on commit 3b09319

Please sign in to comment.