Skip to content

Instantly share code, notes, and snippets.

View coleoguy's full-sized avatar

Heath Blackmon coleoguy

View GitHub Profile
@coleoguy
coleoguy / make.simmap.R
Created September 28, 2016 16:25
modified version of phytools make.simmap
# function creates a stochastic character mapped tree as a modified "phylo" object
# written by Liam Revell 2013, 2014, 2015
# lines 238, 245, 251, and 258 written by Heath Blackmon
make.simmap<-function(tree,x,model="SYM",nsim=1,...){
if(inherits(tree,"multiPhylo")){
ff<-function(yy,x,model,nsim,...){
zz<-make.simmap(yy,x,model,nsim,...)
if(nsim>1) class(zz)<-NULL
return(zz)
}
SimThresh3 <- function(tree, liabilities=F){
## This internal function just determines what the
## observed discrete state is
returnState <- function(x,y){
state <- vector()
# angle A is at origin
# angle B is at 1,0
# angle C is at x,y
# so we can calculate the sides like this
a <- sqrt((x-1)^2 + (y)^2)
library(phytools)
library(evobiR)
tree<-pbtree(n=500, scale=1)
tip.state.full <- SimThresh3(tree, liabilities=T)
ShowTree(tree, tip.state.full[[1]],
cols = c("red", "blue", "green"),
tip.cex=.5, type="fan")
SimThresh3 <- function(tree, liabilities=F){
## This internal function just determines what the
## observed discrete state is
returnState <- function(x,y){
state <- vector()
# angle A is at origin
# angle B is at 1,0
# angle C is at x,y
# so we can calculate the sides like this
a <- sqrt((x-1)^2 + (y)^2)
x <- rnorm(n=100000, mean=1)
mean(x)
x <- sapply(x, StochRound)
mean(x)
@coleoguy
coleoguy / StochRound.R
Created April 30, 2016 19:20
stochastic rounding
StochRound <- function(x){
## extract the decimal portion
q <- abs(x - trunc(x))
## draw a value 0 or 1 with probability
## based on how close we already are
adj <- sample(0:1, size = 1, prob = c(1 - q, q))
## make it negative if x is
if(x < 0) adj <- adj * -1
@coleoguy
coleoguy / blog.parts.R
Created March 9, 2016 22:48
blog.part2
# this time lets set the axis limits to 1:100
plot(0, 0, col = "white", xaxt = "n", yaxt = "n", xlab = "", ylab = "", xlim=c(0,100), ylim=c(0,100))
# now lets try printing a bunch
x <- sample(1:100, size = 100, replace = T)
y <- sample(1:100, size = 100, replace = T)
fl.col <- viridis::viridis(100)
for(i in 1:100){
flower(x = x[i], y = y[i], colx = fl.col[i], cex.x = 1, cex.y = 1.5, cex.f = 1)
}
@coleoguy
coleoguy / example1.R
Created March 9, 2016 22:47
blog.part1
# set up a base plot
plot(0, 0, col = "white", xaxt = "n", yaxt = "n", xlab = "", ylab = "")
flower(x = 0, y = 0, colx = "red", cex.x = .015, cex.y = .035, cex.f = 1)
@coleoguy
coleoguy / flower.R
Last active March 9, 2016 22:45
plots a flower
flower <- function(x, y, colx="red", cex.x=1, cex.y=1, cex.f=1){
petalsx <- cex.x * c(0, 1, 0, -1)
petalsy <- cex.y * c(1, 0, -1, 0)
# this loop prints the petals
for(i in 1:4){
points(x + petalsx[i],
y + petalsy[i],
pch = 16, col = colx,
cex = cex.f)
}
library(shiny)
# Define UI for application
shinyUI(fluidPage(
# Application title
titlePanel("Brownian Motion"),
# Sidebar for user input
sidebarLayout(
sidebarPanel(
# select how many generations to run
sliderInput("gens",