Share your functions & code
- Thread starter trinker
- Start date
- Tags analysis data fundamental graphics
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.
The above might be used given the following rownames in your LaTeX tables (with alt+dragging in notepad++):
The entire latex table can then be re-estimated and replaced in your .tex file with 1 copy+paste command.
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})
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}}} \\
Last edited:
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.
Edit: I added the ability to play multiple notes at the same time and some comments to the code.
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)
}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.
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)
}
}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:
yielding:
Here it is playing with it in R (no knitr):
Here's how you might use it in knitr:
Code:
m1 <- summary(mod)
$p \Sexpr{pform(m1[[3]])}$
Code:
p = .03 #or
p < .01
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)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:
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)
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.
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()@LAzar for some reason it throws up the following error: In readLines(url(y)) : unsupported URL scheme
Haven't really investigated why but want to run it quick to get a better look at the graphical output.
Interesting by the way.
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 connectionInteresting by the way.
