Open cells2numbers opened 7 years ago
Example code to load cell profiler output, track cells and plot windrose plot to visualize migration pattern.
# define column name for track label
strata <- 'TrackObjects_Label'
# read cellprofiler output
profiles <- read_csv("../profiles/201903_testsequence/IdentifyPrimaryObjects.csv") %>%
mutate(Metadata_timepoint = as.numeric(Metadata_timepoint)) %>%
group_by(TrackObjects_Label) %>%
# summarize tracks
tracks <- track(profiles, strata = strata, t_var = "Metadata_timepoint" )
# windrose plot
plot_windrose(tracks)
plot_windrose <- function(population, title_name = "windrose plot", legend = "speed in pixel/frame", dir_res = 30, spd_res = 0.5){
# rotate track angles
population <- tracks %>%
filter(Track_Length > 19) %>%
mutate(Track_Angle = Track_Angle + pi/2) %>% # rotate all angles by 90 degree or pi/2
mutate(Track_Angle = ifelse(Track_Angle > pi, Track_Angle - 2*pi, Track_Angle) ) %>%
mutate(Track_Angle = ifelse(Track_Angle < 0, Track_Angle + 2*pi, Track_Angle) )
h1 <- plot.windrose(
spd = population %>% extract2("Track_Speed"),
dir = (180 * (population$Track_Angle) / pi),
spdmin = 0,
spdmax = max(population$Track_Speed),
spdres = spd_res,
dirres = dir_res,
title_name = title_name,
scale_name = legend
)
}
plot_windrose(tracks)
Code for the complete windrose plot comes from following https://stackoverflow.com/questions/17266780/wind-rose-with-ggplot-r/17266781#17266781
# WindRose.R
#
# following https://stackoverflow.com/questions/17266780/wind-rose-with-ggplot-r/17266781#17266781
require(RColorBrewer)
plot.windrose <- function(data,
spd,
dir,
spdres = 2,
dirres = 30,
spdmin = 2,
spdmax = 20,
spdseq = NULL,
palette = "YlGnBu",
countmax = NA,
debug = 0,
scale_name = "speed in µm/s",
title_name = "Distribution of angle and speed"){
# Look to see what data was passed in to the function
if (is.numeric(spd) & is.numeric(dir)){
# assume that we've been given vectors of the speed and direction vectors
data <- data.frame(spd = spd,
dir = dir)
spd = "spd"
dir = "dir"
} else if (exists("data")){
# Assume that we've been given a data frame, and the name of the speed
# and direction columns. This is the format we want for later use.
}
# Tidy up input data ----
n.in <- NROW(data)
dnu <- (is.na(data[[spd]]) | is.na(data[[dir]]))
data[[spd]][dnu] <- NA
data[[dir]][dnu] <- NA
# figure out the wind speed bins ----
if (missing(spdseq)){
spdseq <- seq(spdmin,spdmax,spdres)
} else {
if (debug >0){
cat("Using custom speed bins \n")
}
}
# get some information about the number of bins, etc.
n.spd.seq <- length(spdseq)
n.colors.in.range <- n.spd.seq - 1
# create the color map
spd.colors <- colorRampPalette(brewer.pal(min(max(3,
n.colors.in.range),
min(9,
n.colors.in.range)),
palette))(n.colors.in.range)
if (max(data[[spd]],na.rm = TRUE) > spdmax){
spd.breaks <- c(spdseq,
max(data[[spd]],na.rm = TRUE))
spd.labels <- c(paste(c(spdseq[1:n.spd.seq-1]),
'-',
c(spdseq[2:n.spd.seq])),
paste(spdmax,
"-",
max(data[[spd]],na.rm = TRUE)))
spd.colors <- c(spd.colors, "grey50")
} else{
spd.breaks <- spdseq
spd.labels <- paste(c(spdseq[1:n.spd.seq-1]),
'-',
c(spdseq[2:n.spd.seq]))
}
data$spd.binned <- cut(x = data[[spd]],
breaks = spd.breaks,
labels = spd.labels,
ordered_result = TRUE)
# clean up the data
data. <- na.omit(data)
# figure out the wind direction bins
dir.breaks <- c(-dirres/2,
seq(dirres/2, 360-dirres/2, by = dirres),
360+dirres/2)
#
dir.labels <- c(paste(360-dirres/2,"-",dirres/2),
paste(seq(dirres/2, 360-3*dirres/2, by = dirres),
"-",
seq(3*dirres/2, 360-dirres/2, by = dirres)),
paste(360-dirres/2,"-",dirres/2))
#dir.labels <- c(paste(0),
# paste(seq(dirres/2, 360-3*dirres/2, by = dirres)),
#"-",
#seq(3*dirres/2, 360-dirres/2, by = dirres)),
# paste(360-dirres/2,"-",dirres/2))
# assign each wind direction to a bin
dir.binned <- cut(data[[dir]],
breaks = dir.breaks,
ordered_result = TRUE)
levels(dir.binned) <- dir.labels
data$dir.binned <- dir.binned
# Run debug if required ----
if (debug>0){
cat(dir.breaks,"\n")
cat(dir.labels,"\n")
cat(levels(dir.binned),"\n")
}
# deal with change in ordering introduced somewhere around version 2.2
if(packageVersion("ggplot2") > "2.2"){
#cat("Hadley broke my code\n")
data$spd.binned = with(data, factor(spd.binned, levels = rev(levels(spd.binned))))
spd.colors = rev(spd.colors)
}
# create the plot ----
p.windrose <- ggplot(data = data,
aes(x = dir.binned,
fill = spd.binned)) +
geom_bar() +
scale_x_discrete(drop = FALSE,
labels = waiver()) +
coord_polar(start = -((dirres/2)/360) * 2*pi) +
scale_fill_manual(name = scale_name,
values = spd.colors,
drop = FALSE) +
theme(axis.title.x = element_blank()) +
labs(title = title_name)
# adjust axes if required
if (!is.na(countmax)){
p.windrose <- p.windrose +
ylim(c(0,countmax))
}
# print the plot
print(p.windrose)
# return the handle to the wind rose
return(p.windrose)
}
Add visualization including