R/plot_helpers.R

Defines functions melt_val_type_data get_avg_stats fdr_line na_sensitive_fun

Documented in fdr_line get_avg_stats melt_val_type_data na_sensitive_fun

# ---------------------------
# plot_helpers.R
# Functions for helping
#	generate_plots.R
# ---------------------------

# If NAs make up >ratio of x, return NA, else return FUN(x, na.rm=TRUE)
#' na_sensitive_fun
#'
#' @export
na_sensitive_fun = function(x, ..., ratio=0.5){
	FUN = ...$FUN

	# Don't just blindly na.rm
	if ((sum(is.na(x)) / length(x)) > ratio){
		return(NA)
	}

	x = x[!is.na(x)]
	return(FUN(x))
}

# FDR
# 1. Order p-values ascending
# 2. Find largest p-value for which p_k <= k/m alpha
# 3. Everything <= p_k
#' fdr_line
#'
#' @export
fdr_line = function(p_vals, alpha=0.05){
	p_vals = p_vals[order(p_vals)]

	# p_k <= (k/length(p_values)) * 0.05
	is_le = (p_vals <= ((1:length(p_vals)/length(p_vals)) * alpha))

	if (!(any(is_le))){
		warning("No p-values are below adjusted cutoff. Returning NA.", call.=TRUE)
		return(NA)
	}

	# grab the largest
	return(max(p_vals[is_le]))
}

# Manually melt
# NOTE: eff, NOT effs
#' get_avg_stats
#'
#' @export
get_avg_stats = function(dists, eff, eff_data, val_name, val_type, do_sd, include_cols, FUN=mean){
	if (length(eff) > 1){
		warning("Please pass a single efficacy value. Aborting...", .call=TRUE)
		quit()
	}

	avg_data = data.frame()
	for (d in 1:length(dists)){
		to_add = aggregate(as.formula(paste(val_name, '~', val_type, sep=' ')),
				   eff_data[eff_data[,"dist"] == dists[d],],
				   FUN=na_sensitive_fun, ...=list(FUN=FUN))

		if (do_sd){
			add_sd = aggregate(as.formula(paste(val_name, '~', val_type, sep=' ')),
					   eff_data[eff_data[,"dist"] == dists[d],],
					   FUN=na_sensitive_fun, ...=list(FUN=sd))
			colnames(add_sd) = c(val_type, "avg_sd")

			to_add = cbind(to_add, add_sd)
			to_add = to_add[,!duplicated(colnames(to_add))]
		}

		if (include_cols){
			to_add = cbind(dist=dists[d], eff=eff, to_add)
		}

		# NOTE: jank way of detecting volcano plot or not
		if (!do_sd){
			colnames(to_add)[colnames(to_add) == val_name] = paste("avg_", val_name, sep='')

		} else {
			colnames(to_add)[colnames(to_add) == val_name] = "avg_value"
		}

		avg_data = rbind(avg_data, to_add)
	}

	return(avg_data)
}

#' melt_val_type_data
#'
#' @export
melt_val_type_data = function(x, val_name, val_type){
	to_keep = grepl(val_name, colnames(x)) |
		colnames(x) %in% c("eff", "dist")

	tmp = x[,to_keep]
	tmp_cols = colnames(tmp)[grepl(val_name, colnames(tmp))]

	# Convert wide to long
	tmp = melt(tmp, id.vars=c("dist", "eff"), measure.vars=tmp_cols,
			variable_name=val_type)
	colnames(tmp)[colnames(tmp) == "value"] = val_name

	return(tmp)
}
kmorrisongr/ksmthesis documentation built on Oct. 5, 2020, 6:41 a.m.