Skip to content

Instantly share code, notes, and snippets.

@jwaage
jwaage / sinecurve.R
Created September 15, 2017 08:10
Sine curve function and plot
# Season test
# Do a sine function with arbitrary phase shift and plot the results with confidence band.
library(tidyverse)
library(broom)
set.seed(1)
x <- runif(200, 1, 365) %>% round(0)
y <- sin(2*pi*x/365) + rnorm(200, 0, 0.5)
yshift <- sin(2*pi*x/365 - 0.3*pi) + rnorm(200, 0, 0.5)
@jwaage
jwaage / adonis_table.R
Created August 18, 2017 14:14
Adonis Table
adonis_table <- function(phy, d, ...){
subs <- function(dist, index) as.dist(as.matrix(dist)[index,index])
contcheck <- function(var){
suppressWarnings(
if(all(is.na(as.numeric(var))) | length(unique(var))<10){
return(var)
} else {
return(as.numeric(var))
}
)
@jwaage
jwaage / string.R
Created October 12, 2016 14:06
StringDB PPI grabber
# Libraries
# ------------------------------------------------------------------------------
library("RCurl")
library("bitops")
# Define functions
# ------------------------------------------------------------------------------
getInteractions = function(id,score,limit,species=9606,print=FALSE){
@jwaage
jwaage / ggphylobarplot.R
Created June 28, 2016 12:40
Phyloseq Barplot, customizable ggplot2
ggphylobarplot <- function(physobj,
barlevel = "lowest",
colorlevel = "phylum",
splitvar = "Time",
splitlevelorder = c("1w", "1m", "3m"),
splitlevellabels = c("1 week", "1 month", "3 months"),
sortlevel = "1w",
mracut = .009,
agglomerate = FALSE,
agglomeratelevel = "genus",
@jwaage
jwaage / CV_splsda.R
Created June 16, 2016 21:26
mixOmics splsda - CV function for component and variable selection
getBestFit <- function(x, y, ..., ncomp = 5, nfolds = 5, minvar = 1, maxvar = ncol(x), direction = ">", seed = 42, run = NULL){
set.seed(seed)
auc <- expand.grid(AUC = 0, ncomp = 1:ncomp, nvar = minvar:maxvar, stringsAsFactors = FALSE)
folds <- sample.int(nfolds, size = nrow(x), replace = TRUE)
for(i in minvar:maxvar){
pred <- matrix(rep(NA,nrow(x)*ncomp), ncol = ncomp)
for(fold in 1:nfolds){
train <- folds != fold
test <- folds == fold
fit <- splsda(x[train,], y[train], keepX = rep(i,ncomp), ncomp = ncomp)
@jwaage
jwaage / CV_contrib.R
Created June 16, 2016 21:22
mixOmics splsda - examine variable contributions across CV folds
getBugsCV <- function(x, y, ncomp = 5, nfolds = 5, nvar = 5, direction = ">", seed = 42, nreps = 1){
set.seed(seed)
sigbugs <- list(nfolds*nreps)
for(rep in 1:nreps){
folds <- sample.int(nfolds, size = nrow(x), replace = TRUE)
for(fold in 1:nfolds){
train <- folds != fold
test <- folds == fold
fit <- splsda(x[train,], y[train], keepX = rep(nvar,ncomp), ncomp = ncomp)
sigbugs[[(rep-1)*nfolds + fold]] <- names(fit$loadings$X[,1])[fit$loadings$X[,1] != 0]
@jwaage
jwaage / lowest.R
Last active June 17, 2016 22:11
Create "lowest" tax rank for phyloseq
lowest <- function(tax){
returnPhyseq <- FALSE
if(class(tax) == "phyloseq"){
physeq <- tax
tax <- tax_table(physeq)
returnPhyseq <- TRUE
}
tax <- tax[, colnames(tax) != "lowest"]
tax <- tax %>% as("matrix")
low <- tax[,ncol(tax)] %>% as.vector
@jwaage
jwaage / ggroc.R
Last active March 5, 2021 20:34
ROC curve using ggplot2 and pROC
ggroc <- function(roc, showAUC = TRUE, interval = 0.2, breaks = seq(0, 1, interval)){
require(pROC)
if(class(roc) != "roc")
simpleError("Please provide roc object from pROC package")
plotx <- rev(roc$specificities)
ploty <- rev(roc$sensitivities)
ggplot(NULL, aes(x = plotx, y = ploty)) +
geom_segment(aes(x = 0, y = 1, xend = 1,yend = 0), alpha = 0.5) +
geom_step() +
@jwaage
jwaage / get_loadings.R
Created March 3, 2016 14:08
Phyloseq PCoA Loadings
get_loadings <- function(physeq, ord, method = c("pearson", "spearman"), scaled=TRUE, cor = FALSE, axes=1:2, taxranks=6:7, top = NULL){
if(!all(rownames(ord$vectors) == sample_names(physeq)))
stop("Taxa names do not match")
counts <- t(as(otu_table(physeq),"matrix"))
scores <- ord$vectors
cov <- cov(counts, scores, method=method[1])
if(cor){
cov <- cor(counts, scores, method=method[1])
}
if(scaled){
@jwaage
jwaage / googlesheet.R
Created February 3, 2016 12:55
Get data from google sheet
library("googlesheets")
# Forbind til google. Efter denne kommando åbner browseren, og du skal logge ind på google.
gs_ls()
# Herefter registerer du det pågældende ark med enten ark-navnet, eller som her, ark-nøglen (står i URLen)
sheet_key <- "15WDkfscZ-3xuSE33WRmj5PIpC-q5Eto8wUlWqlv-x2w"
sheet_object <- gs_key(sheet_key)
sheet_csv <- gs_read_csv(sheet_object, "Fanenavn")