Share your functions & code

Dason

Ambassador to the humans
Does alarm() work for you though? And I'll admit I only tested that function on my laptop and since I'm no good at Windows command line stuff and I basically just stole that call to system from an example page... I never had high hopes for it being able to work on any Windows system.
 

trinker

ggplot2orBust
Sad :( no but maybe I turned it off. My compter came complete with bells and whistles literly. It annoyed me and all noises soon were turned off. YEah I wanted a trumpet to sound every time I opened a web page.
 
In my casual RA work I have to make 10s of results tables each month and then put them into LaTeX. Often I will have to re-estimate them many times if I want to change a tiny thing. So I made a function which is to be used in conjunction with xtable and notepad++ (for alt+dragging --> copy+paste into WinEdt or equivalent). not useful for most of you probably.

Code:
fun_breakup <- function(matrix,blocksizes,gapsizes)
{
require(miscTools)
counter <- 0
for(i in 1:length(blocksizes))
for(j in 1:gapsizes[i])
{{
	counter <- counter + 1
	if(j==1){ matrix <- insertRow(matrix,sum(blocksizes[1:i])+counter) }
    if(j>1) { matrix <- insertRow(matrix,sum(blocksizes[1:i])+counter) }	
}}

return(matrix)
}
Code:
matrix <- cbind(rnorm(50),rnorm(50))
blocks <- c(5,5,3,5,3,5,3,5,3,5,3)   
gaps <- c(2,1,2,1,2,1,2,1,2,1,2)
require(xtable)
print(xtable(  fun_breakup(matrix,blocks,gaps)  ),type="latex",sanitize.text.function=function(x){x})
The above might be used given the following rownames in your LaTeX tables (with alt+dragging in notepad++):

Code:
\multicolumn{5}{@{}l}{\emph{\underline{Variable 1}}} \\
$\tau = 0.05$  .......
$\tau = 0.10$  .......
$\tau = 0.50$  .......
$\tau = 0.90$  .......
$\tau = 0.95$  .......
\\
\multicolumn{5}{@{}l}{\emph{\underline{Variable 2}}} \\
$\tau = 0.05$ .......
$\tau = 0.10$ .......
$\tau = 0.50$ .......
$\tau = 0.90$ .......
$\tau = 0.95$ .......
\\
$F:\tau \in \{0.05,0.50\}$.........
$F:\tau \in \{0.5,0.95\}$ .........
$F:\tau \in \{0.05,0.95\}$.........
\\
\multicolumn{5}{@{}l}{\emph{\underline{Variable 3}}} \\
The entire latex table can then be re-estimated and replaced in your .tex file with 1 copy+paste command.
 
Last edited:

Dason

Ambassador to the humans
After watching this today I decided to see if I could make some music in R. Turns out you can (sort of). R can't actually play the resulting file but you can create .wav files easily enough. See if you can guess the song before actually running the code. You'll need the tuneR and XML packages installed. I also tried to make it both Windows and Linux friendly but was too lazy to make it work for macs. I think all you would need to do is change "xdg-open" to "open" if you're on a mac though. Enjoy.

Code:
library(tuneR)
library(XML)

# Grab note frequencies off the web
notes <- readHTMLTable("http://www.phy.mtu.edu/~suits/notefreqs.html")[[1]]
freq <- as.numeric(as.character(notes[,2]))
# Make life easier by only referring to sharps
names(freq) <- substr(notes[,1], 1, 3)

# Only takes two values currently and they are assummed
# to be the same length.  Constructs a chord out of the
# two Waves
chord <- function(x, y){
    out <- x
    out@left <- x@left + y@left
    out <- normalize(out, unit = as.character(x@bit))
    return(out)
}

# Only takes a single character string input and
# returns a Wave corresponding to that note
# convert_to_sine("A4")
convert_to_sine <- function(note, duration = 1/3, wave = sine){
    if(note == "rest"){
        return(silence(duration = duration, bit = 16, xunit = "time"))
    }
    # Check for chords - These have periods in them
    if(grepl(".", note, fixed = TRUE)){
        notes <- strsplit(note, "\\.")[[1]]
        sine_waves <- lapply(notes, convert_to_sine)
        return(chord(sine_waves[[1]], sine_waves[[2]]))
    }
    # If just a single note...
    wave(freq[note], bit = 16, xunit = "time", duration = duration)
}

# Song Definition
intro  <- c("E6", "D#6", 
            "E6", "D#6", "E6",  "B5",  "D6",  "C6")
v1     <- c("A3.A5", "E4.A5",  "A4",  "C5",  "E5",  "A5",
            "E3.B5", "E4.B5",  "G#4", "E5",  "G#5", "B5",
            "A3.C6", "E4.C6",  "A4",  "E5",  "E6",  "D#6",
            "E6", "D#6", "E6",  "B5",  "D6",  "C6",
            "A3.A5", "E4.A5",  "A4",  "C5",  "E5",  "A5",
            "E3.B5", "E4.B5",  "G#4", "E5",  "C6",  "B5")
v1e1   <- c("A3.A5", "E4.A5",  "A4",  "rest", intro)
v1e2   <- c("A3.A5", "E4.A5",  "A4",  "B5",  "C6",  "D6")
v2     <- c("C4.E6", "G4.E6",  "C5.E6",  "G5",  "F6",  "E6",
            "G4.D6", "G4.D6",  "B4.D6",  "F5",  "E6",  "D6",
            "A3.C6", "E4.C6",  "A4.C6",  "E5",  "D6",  "C6",
            "E3.B5", "E4.B5",  "E5",  "E5",  "E6",  "E5",
            "E6", "E5",  "E6",  "D#6", "E6",  "D#6",
            "E6", "D#6", "E6",  "D#6", "E6",  "D#6",
            "E6", "D#6", "E6",  "B5",  "D6",  "C6",
            v1)
v2e1    <- c("A3.A5", "E4.A5",  "A4",  "B5",  "C6",  "D6")
v2e2    <- c("A3.A5", "A3.A5")
song_notes <- c(intro, v1, v1e1, v1, v1e2, v2, v2e1, v2, v2e2)

# Turn the definition into Wave objects
sine_waves <- lapply(song_notes, convert_to_sine)
# Smooth out the transitions
smoothed_waves <- suppressWarnings(lapply(sine_waves, prepComb))
# Bind the Waves together
song <- do.call("bind", smoothed_waves)

# Actually play the song
if(.Platform['OS.type'] == "windows"){
    play(song)
}else{
    tmpfile <- if(exists("tmpfile")){tmpfile}else{tempfile(fileext = ".wav")}
    writeWave(song, file = tmpfile)
    com <- paste("xdg-open", tmpfile)
    system(com)
}
Edit: I added the ability to play multiple notes at the same time and some comments to the code.
 

Dason

Ambassador to the humans
Ok time for one more song...

I've also modified the chord function because I needed to be able to take up to 4 notes at a time for this song. I renamed a few things because the new names seem more appropriate.
Code:
library(tuneR)
library(XML)

# Grab note frequencies off the web
notes <- readHTMLTable("http://www.phy.mtu.edu/~suits/notefreqs.html")[[1]]
freq <- as.numeric(as.character(notes[,2]))
# Make life easier by only referring to sharps
names(freq) <- substr(notes[,1], 1, 3)

# The input is multiple Waves that are combined together
chord <- function(...){
    out <- list(...)[[1]]
    out@left <- rowSums(sapply(list(...), slot, name = "left"))
    out <- normalize(out, unit = as.character(out@bit))
    return(out)
}

# Only takes a single character string input and
# returns a Wave corresponding to that note
# convert_to_sine("A4") # Gives a single A4 note
# convert_to_sine("A4.E5") # Gives a chord combining A4 and E5
# convert_to_sine("2_A4.E5") # Gives the above chord with twice the duration
convert_to_wave <- function(note, duration = 1/4, wave = sine){
    # Check for a modified duration
    if(grepl("_", note, fixed = TRUE)){
        tmp <- strsplit(note, "_")[[1]]
        duration_modifier <- as.numeric(tmp[1])
        return(convert_to_wave(tmp[2], duration = duration*duration_modifier, wave = wave))
    }
    
    if(note == "rest"){
        return(silence(duration = duration, bit = 16, xunit = "time"))
    }
    
    # Check for chords - These have periods in them
    if(grepl(".", note, fixed = TRUE)){
        notes <- strsplit(note, "\\.")[[1]]
        waves <- lapply(notes, convert_to_wave, duration = duration, wave = wave)
        return(do.call(chord, waves))
    }
    # If just a single note...
    wave(freq[note], bit = 16, xunit = "time", duration = duration)
}

# Song Definition
base_duration <- 1/4
intro <- c("C4", "C4.G4", "C4.E4", "C4.G4", "C4", "C4.G4", "C4.E4", "C4.G4",
           "C4", "C4.G4", "C4.F4", "C4.G4", "C4", "C4.G4", "C4.F4", "C4.G4",
           "C4", "C4.G4", "C4.E4", "C4.G4", "C4", "C4.G4", "C4.E4", "C4.G4",
           "C4", "C4.G4", "C4.F4", "C4.G4", "C4", "C4.G4", "C4.F4", "C4.G4")
m1    <- c("C4.E5", "G4.E5", "E4.E5", "G4.E5", "C4.F5", "G4.F5", "E4.G5", "G4.G5",
           "B3.G5", "G4.G5", "E4.F5", "G4.F5", "B3.E5", "G4.E5", "D4.D5", "G4.D5",
           "A3.C5", "G4.C5", "E4.C5", "G4.C5", "A3.D5", "G4.D5", "E4.E5", "G4.E5",
           "G3.E5", "G4.E5", "D4.E5", "D4.D5", "G3.D5", "G4.D5", "D4.D5", "G4.D5")
m2    <- c("C4.E5", "G4.E5", "E4.E5", "G4.E5", "C4.F5", "G4.F5", "E4.G5", "G4.G5",
           "B3.G5", "G4.G5", "E4.F5", "G4.F5", "B3.E5", "G4.E5", "D4.D5", "G4.D5",
           "A3.C5", "G4.C5", "E4.C5", "G4.C5", "A3.D5", "G4.D5", "E4.E5", "G4.E5",
           "G3.D5", "G4.D5", "D4.D5", "D4.C5", "C4.C5", "G4.C5", "E4.C5", "G4.C5")
m3    <- c("B3.D5", "G4.D5", "D4.D5", "G4.D5", "C4.E5", "G4.E5", "E4.C5", "G4.C5",
           "B3.D5", "G4.D5", "G4.E5", "G4.F5", "C4.E5", "G4.E5", "E4.C5", "G4.C5",
           "B3.D5", "G4.D5", "G4.E5", "G4.F5", "E3.E5", "G#4.E5", "G#4.D5", "G#4.D5",
           "2_A3.E4.A5.C5", "2_F#3.D4.A4.D5", "4_G3.D4.G4")
m4    <- c("C4.E5", "G4.E5", "E4.E5", "G4.E5", "C4.F5", "G4.F5", "E4.G5", "G4.G5",
           "B3.G5", "G4.G5", "E4.F5", "G4.F5", "B3.E5", "G4.E5", "D4.D5", "G4.D5",
           "A3.C5", "G4.C5", "E4.C5", "G4.C5", "A3.D5", "G4.D5", "E4.E5", "G4.E5",
           "G3.D5", "G4.D5", "D4.D5", "D4.C5", "4_C4.E4.G4.C5")

song_notes <- c(intro, m1, m2, m3, m4)
song_notes <- m3

# Turn the definition into Wave objects
waves <- lapply(song_notes, convert_to_wave, duration = base_duration)
# Smooth out the transitions
smoothed_waves <- suppressWarnings(lapply(waves, prepComb))
# Bind the Waves together
song <- do.call("bind", smoothed_waves)

# Actually play the song
if(.Platform['OS.type'] == "windows"){
    play(song)
}else{
    tmpfile <- if(exists("tmpfile")){tmpfile}else{tempfile(fileext = ".wav")}
    writeWave(song, file = tmpfile)
    com <- paste("xdg-open", tmpfile)
    system(com)
}
 
Function for turning a daily time series to a monthly time series. "dates" is a vector of daily character strings "1/01/1988", "2/01/1988",... for example. Set over="first" to make the monthly time series from the first day of each month, and over="last" to make the time series to the last day of each month.

Code:
daily_to_monthly <- function(over=c("first","last"),dates)
{
	splits <- strsplit(dates,"/")
	splits <- matrix(do.call(rbind,splits),ncol=3)
	
	last <- sapply(2:nrow(splits), function(i) { as.numeric(splits[i,2])!=as.numeric(splits[i-1,2]) })
	if(as.numeric(splits[1,1])<3){
		first <- c(TRUE,last) 
	} else {
		first <- c(FALSE,last)
	}

if(over=="first"){
	return(first)
} else {
	return(last)
}	
}
 

trinker

ggplot2orBust
I got a little one but it's becoming handy. When you use knitr or sweave and want to include pvalues in the code some may be so close to one that if you round 2 digits you get p = 0 which bothers me. I decided to write a function that takes the pvalue as an argument and optionally a cutoff value (default is .01) and a digits argument (default is null and matches digits found in the cutoff). it automatically switches to < if the p value is less than the cutoff. simple but was kinda nice in the knitr assignment I just did.

Here's how you might use it in knitr:

Code:
m1 <- summary(mod)
$p \Sexpr{pform(m1[[3]])}$
yielding:

Code:
p = .03 #or
p < .01
Here it is playing with it in R (no knitr):

Code:
pform <- function(x, cutoff = .01, digits = NULL) {
    z <- strsplit(as.character(cutoff), "\\.")[[c(1, 2)]]
    w <- paste0(".", z)
    if (x < cutoff) {
        paste0("< ", w)
    } else {
        if (is.null(digits)) digits <- nchar(z)
        paste0("= .", strsplit(as.character(round(x, digits)), "\\.")[[c(1, 2)]])
    }
}#end of pform function

#EXAMPLE:
pform(.001)
pform(.3)
pform(.0013, cutoff = .00001, digits = 3)
 

trinker

ggplot2orBust
Try it and you'll see it explicitly takes care of the = or < sign for you. Just a slight time saver butt if you repeat it enough...

I don't think format.pval does that does it? If so I've done it again; recreated what already was :)
 

trinker

ggplot2orBust
Linux users this isn't for you sorry :(

If you use dput but are lazy like me here's a function to copy the output from dput to your clipboard if you are a mac or Win machine:

Code:
repex <- function(x) {
    v <- dput(x)
    z <- capture.output(dput(x))
    if (Sys.info()["sysname"] == "Windows") {
        writeClipboard(z, format = 1)
    }
    if (Sys.info()["sysname"] == "Darwin") {           
        j <- pipe("pbcopy", "w")                       
        writeLines(z, con = j)                               
        close(j)                                    
    }     
}

repex(mtcars)
 

Lazar

Phineas Packard
Someone on facebook claimed that suicide rates were much much higher in the west. I thought BS so did some webscraping to see. Here is the resulting graph:

And here is the code. I am sure others would find a better way to match more countries and make the code neater but this was fun.
Code:
#Suicide data from world health organisaion
library(XML)
x <- 'http://www.who.int/mental_health/prevention/suicide_rates/en/'
x<- readHTMLTable(x)
#The table consists a bunch of empty lines so I got rid of them.
suicide <- cbind(x[[1]][1:105,1:4])
#Get rid of white space
suicide[,1]<- gsub('[[:space:]]', '', suicide[,1])

#CIA data for GDP per capita
y<-'https://www.cia.gov/library/publications/the-world-factbook/rankorder/rawdata_2004.txt'
y <- readLines(url(y))
#Parse data to place in dataframe
y <- gsub('\\t','', y)
y <- gsub('^[0-9]+','', y)
gdp <- strsplit(y, '\\$')
gdp<- data.frame(do.call(rbind, gdp))
gdp[,2] <- as.numeric(gsub(',','', gdp[,2]))
names(gdp)<- c("Country", 'GPD')
#Change countries to upper case to match with suicide rates
gdp[,1]<- toupper(gdp[,1])

#Merge data
data <- merge(gdp, suicide, by ='Country')
#Change relavent variables from character to numeric
data[,4]<- as.numeric(data[,4])
data[,5]<- as.numeric(data[,5])
#Mean suicide rates across males and females
#Based on assumption that genders are roughly evenly split by country
data$total <- apply(data[,4:5],1, mean)
#Correlations
cor(data[,c(2,4:6)])
#For graphing need a new dataset
Male = data[,1:4]
Male$Gender <- rep('blue', nrow(Male))
names(Male)[4]<- 'suicide'
Female = data[,c(1:3,5)]
Female$Gender <- rep('pink', nrow(Female))
names(Female)[4]<- 'suicide'
dataVis <- rbind(Male, Female)

#svg('C:/Users/Philip/ubuntu1/facebookGraph.svg', width = 11, height = 7)
plot(dataVis[,2],dataVis$suicide, col = as.character(dataVis$Gender),
     pch = '', ylab='Suicides per 10,000', xlab='GDP per capita',
     main='Relationship Between National Wealth and Suicide \n By Gender')
abline(lm(dataVis[dataVis$Gender=='blue',4]~dataVis[dataVis$Gender=='blue',2]),
       col='blue')
abline(lm(dataVis[dataVis$Gender=='pink',4]~dataVis[dataVis$Gender=='pink',2]),
       col='pink')
text(dataVis[,2],dataVis$suicide,dataVis$Country, 
     col = as.character(dataVis$Gender), cex=.5)
legend(70000, 80, c('Male','Female'), text.col=c('blue', 'pink'), cex=.8,
       box.lwd = 0, ,box.col = "white",bg = "white")
#dev.off()
 

trinker

ggplot2orBust
@LAzar for some reason it throws up the following error: In readLines(url(y)) : unsupported URL scheme
Code:
> y<-'http://www.cia.gov/library/publications/the-world-factbook/rankorder/rawdata_2004.txt'
> y <- readLines(url(y))
Error in readLines(url(y)) : cannot open the connection
Haven't really investigated why but want to run it quick to get a better look at the graphical output.

Interesting by the way.
 

Dason

Ambassador to the humans
trinker it looks like you changed it from https to http. The https is whats causing the problem - R can't directly read from https using the methods Lazar used so I'm interested in how he got that to work.
 

trinker

ggplot2orBust
Yeah I assumed at first it was the https so I made the switch but still get the same error. I'm curious on two accounts now:

1. how to get it to work
2. what the produced visual will tell
 

Dason

Ambassador to the humans
One quick and nonbrainy solution. Download it in your system and give the change to y="\home\antibot\rawdata_2004.txt"
That's what I ended up doing because I wanted to see it plotted on a log-log scale but I'm still interested in how Lazar got it to work.