Skip to content

Commit ce5ad6a

Browse files
author
jaywarrick
committed
initial commit
1 parent c3e53d3 commit ce5ad6a

6 files changed

+593
-0
lines changed

‎drawLogicleAxis.R

+22
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,22 @@
1+
drawLogicleAxis <- function(ticks=c(-10,1,10,100,1000,10000,100000,1000000), tickLabels=ticks, axisNum=1, transition=1, linLogRatio=1)
2+
{
3+
#tickLabels <- tickLabels[1:length(ticks)]
4+
# Create the tick labels and tick locations
5+
# for(i in 1:length(ticks))
6+
# {
7+
# if(ticks[i] == transition)
8+
# {
9+
# tickLabels[i] = paste(as.character(tickLabels[i]), '*', sep ='')
10+
# }
11+
# }
12+
ticks <- logicle(ticks, transition, linLogRatio)
13+
if(axisNum == 2)
14+
{
15+
axis(axisNum, at=ticks, labels=tickLabels, las=2)
16+
}
17+
else
18+
{
19+
axis(axisNum, at=ticks, labels=tickLabels)
20+
}
21+
22+
}

‎getDoubleStats.R

+48
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,48 @@
1+
getDoubleStats <- function(x=1:10, y=1:10, xThresh=5, xCross=0, yThresh=5, yCross=0)
2+
{
3+
ret <- getSingleStats(TRUE, x, y, xThresh, xCross);
4+
ret$values <- getSingleStats(FALSE, x, y, yThresh, yCross);
5+
6+
total <- ret$'n'
7+
8+
# Define subpopulations
9+
plusplus <- ((x > (y*yCross + xThresh)) & (y > (x*xCross + yThresh)))
10+
minusminus <- ((x <= (y*yCross + xThresh)) & (y <= (x*xCross + yThresh)))
11+
plusminus <- ((x > (y*yCross + xThresh)) & (y <= (x*xCross + yThresh)))
12+
minusplus <- ((x <= (y*yCross + xThresh)) & (y > (x*xCross + yThresh)))
13+
14+
# Calculate stats for each subpopulation
15+
# ++
16+
ret$'XY % ++' <- 100 * (sum(plusplus) / (total))
17+
ret$'XY n ++' <- sum(plusplus)
18+
ret$'XY meanX ++' <- mean(x[plusplus])
19+
ret$'XY meanY ++' <- mean(y[plusplus])
20+
ret$'XY sdX ++' <- sd(x[plusplus])
21+
ret$'XY sdY ++' <- sd(y[plusplus])
22+
23+
# --
24+
ret$'XY % --' <- 100 * (sum(minusminus) / (total))
25+
ret$'XY n --' <- sum(minusminus)
26+
ret$'XY meanX --' <- mean(x[minusminus])
27+
ret$'XY meanY --' <- mean(y[minusminus])
28+
ret$'XY sdX --' <- sd(x[minusminus])
29+
ret$'XY sdY --' <- sd(y[minusminus])
30+
31+
# -+
32+
ret$'XY % -+' <- 100 * (sum(minusplus) / (total))
33+
ret$'XY n -+' <- sum(minusplus)
34+
ret$'XY meanX -+' <- mean(x[minusplus])
35+
ret$'XY meanY -+' <- mean(y[minusplus])
36+
ret$'XY sdX -+' <- sd(x[minusplus])
37+
ret$'XY sdY -+' <- sd(y[minusplus])
38+
39+
# +-
40+
ret$'XY % +-' <- 100 * (sum(plusminus) / (total))
41+
ret$'XY n +-' <- sum(plusminus)
42+
ret$'XY meanX +-' <- mean(x[plusminus])
43+
ret$'XY meanY +-' <- mean(y[plusminus])
44+
ret$'XY sdX +-' <- sd(x[plusminus])
45+
ret$'XY sdY +-' <- sd(y[plusminus])
46+
47+
return(ret);
48+
}

‎getSingleStats.R

+55
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,55 @@
1+
getSingleStats <- function(isVerticalThreshold, x, y, cross, thresh)
2+
{
3+
ret <- hash()
4+
ret$n <- length(x)
5+
ret$meanX <- mean(x)
6+
ret$meanY <- mean(y)
7+
ret$sdX <- sd(x)
8+
ret$sdY <- sd(y)
9+
10+
total <- ret$n
11+
12+
if(isVerticalThreshold)
13+
{
14+
# Define subpopulations
15+
plus <- (x > (y*cross + thresh))
16+
minus <- (x <= (y*cross + thresh))
17+
18+
# Calculate stats for each subpopulation
19+
ret$'X % +' <- 100 * (sum(plus) / (total));
20+
ret$'X n +' <- sum(plus)
21+
ret$'X meanX +' <- mean(x[plus])
22+
ret$'X sdX +' <- sd(x[plus])
23+
ret$'X meanY +' <- mean(y[plus])
24+
ret$'X sdY +' <- sd(y[plus])
25+
26+
ret$'X % -' <- 100 * (sum(minus) / (total))
27+
ret$'X n -' <- sum(minus)
28+
ret$'X meanX -' <- mean(x[minus])
29+
ret$'X sdX -' <- sd(x[minus])
30+
ret$'X meanY -' <- mean(y[minus])
31+
ret$'X sdY -' <- sd(y[minus])
32+
}
33+
else
34+
{
35+
# Define subpopulations
36+
plus <- (y > (x*cross + thresh))
37+
minus <- (y <= (x*cross + thresh))
38+
39+
# Calculate stats for each subpopulation
40+
ret$'Y % +' <- 100 * (sum(plus) / (total))
41+
ret$'Y n +' <- sum(plus)
42+
ret$'Y meanX +' <- mean(x[plus])
43+
ret$'Y sdX +' <- sd(x[plus])
44+
ret$'Y meanY +' <- mean(y[plus])
45+
ret$'Y sdY +' <- sd(y[plus])
46+
47+
ret$'Y % -' <- 100 * (sum(minus) / (total))
48+
ret$'Y n -' <- sum(minus)
49+
ret$'Y meanX -' <- mean(x[minus])
50+
ret$'Y sdX -' <- sd(x[minus])
51+
ret$'Y meanY -' <- mean(y[minus])
52+
ret$'Y sdY -' <- sd(y[minus])
53+
}
54+
return(ret)
55+
}

‎logicle.R

+10
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,10 @@
1+
logicle <- function(x, transitionPoint, linearUnitsPerOrder)
2+
{
3+
linearDifference = x - transitionPoint;
4+
valsToAdjust <- linearDifference <= 0;
5+
ordersDifferenceOnDisplay = linearDifference / linearUnitsPerOrder;
6+
linearDifference[valsToAdjust] <- transitionPoint * 10^(ordersDifferenceOnDisplay[valsToAdjust]);
7+
x[valsToAdjust] <- transitionPoint*10^(-1*((transitionPoint-x[valsToAdjust])/linearUnitsPerOrder))
8+
9+
return(x)
10+
}

‎plotFACSPercentages.R

+164
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,164 @@
1+
# rm(list = ls())
2+
3+
# xmin=-100
4+
# xmax=500000
5+
# ymin=-100
6+
# ymax=500000
7+
# xThresh=74.79895
8+
# yThresh=168.0746
9+
# xCross=0
10+
# yCross=0.002
11+
# xTransition=xThresh
12+
# yTransition=yThresh
13+
# xLinLogRatio=200
14+
# yLinLogRatio=xLinLogRatio
15+
# xTicks=c(0,100,1000,10000,100000,1000000)
16+
# yTicks=xTicks
17+
18+
source('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R FACS Plotting Functions/logicle.R')
19+
source('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R FACS Plotting Functions/plotFACS.R')
20+
source('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R FACS Plotting Functions/getSingleStats.R')
21+
source('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R FACS Plotting Functions/getDoubleStats.R')
22+
source('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R FACS Plotting Functions/drawStats.R')
23+
source('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R FACS Plotting Functions/drawLogicleAxis.R')
24+
source('/Users/jaywarrick/GoogleDrive/SingleCell/Figures/plotHelperFunctions.R')
25+
library(flowCore)
26+
library(flowViz)
27+
library(MASS)
28+
library(RColorBrewer)
29+
library(foreign)
30+
library(hash)
31+
32+
FACSFolder <- '/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/Compiled Data Folders/FCS Data/Compiled/'
33+
setwd('/Users/jaywarrick/GoogleDrive/SingleCell/FACS Data/R Plots/')
34+
35+
36+
thresh <- list()
37+
med <- list()
38+
files <- c('WT_6_B.fcs', 'WT_12_C.fcs', 'WT_S_A.fcs', 'M51R_6_A.fcs', 'M51R_12_C.fcs', 'M51R_S_B.fcs', 'Mock_6_C.fcs', 'Mock_12_C.fcs', 'Mock_S_B.fcs')
39+
data <- read.FCS(paste(FACSFolder, 'Mock_6_C.fcs', sep=''), transformation=FALSE)
40+
temp <- data.frame(G=exprs(data$'488 B 530/30-A')[,1], R=exprs(data$'561 D 610/20-A')[,1])
41+
med$Early <- c(median(temp$G), median(temp$R))
42+
mad <- c(mad(temp$G), mad(temp$R))
43+
thresh$Early <- 3*mad
44+
45+
data <- read.FCS(paste(FACSFolder, 'Mock_12_C.fcs', sep=''), transformation=FALSE)
46+
temp <- data.frame(G=exprs(data$'488 B 530/30-A')[,1], R=exprs(data$'561 D 610/20-A')[,1])
47+
med$Middle <- c(median(temp$G), median(temp$R))
48+
mad <- c(mad(temp$G), mad(temp$R))
49+
thresh$Middle <- 3*mad
50+
51+
data <- read.FCS(paste(FACSFolder, 'Mock_S_B.fcs', sep=''), transformation=FALSE)
52+
temp <- data.frame(G=exprs(data$'488 B 530/30-A')[,1], R=exprs(data$'561 D 610/20-A')[,1])
53+
med$Late <- c(median(temp$G), median(temp$R))
54+
mad <- c(mad(temp$G), mad(temp$R))
55+
thresh$Late <- 3*mad
56+
57+
times2 <- c(6,12,18)
58+
results1 <- list()
59+
for(f in files)
60+
{
61+
data <- read.FCS(paste(FACSFolder, f, sep=''), transformation=FALSE)
62+
temp <- data.frame(G=exprs(data$'488 B 530/30-A')[,1], R=exprs(data$'561 D 610/20-A')[,1])
63+
if(length(grep('6', f, fixed=TRUE) > 0))
64+
{
65+
gateX <- thresh$Early[1]
66+
gateY <- thresh$Early[2]
67+
temp$G <- temp$G - med$Early[1]
68+
temp$R <- temp$R - med$Early[2]
69+
}
70+
if(length(grep('12', f, fixed=TRUE) > 0))
71+
{
72+
gateX <- thresh$Middle[1]
73+
gateY <- thresh$Middle[2]
74+
temp$G <- temp$G - med$Middle[1]
75+
temp$R <- temp$R - med$Middle[2]
76+
}
77+
if(length(grep('S', f, fixed=TRUE) > 0))
78+
{
79+
gateX <- thresh$Late[1]
80+
gateY <- thresh$Late[2]
81+
temp$G <- temp$G - med$Late[1]
82+
temp$R <- temp$R - med$Late[2]
83+
}
84+
stats <- getDoubleStats(temp$G, temp$R, xThresh=gateX, xCross=0, yThresh=gateY, yCross=0.002)
85+
results1[[f]] <- stats
86+
}
87+
88+
89+
######## Plot 1 ###############
90+
fun <- function(hashObj, hashKey)
91+
{
92+
return(hashObj[[hashKey]])
93+
}
94+
pdf('WT_FACSPercentages_Plot.pdf',width=1.1*W,height=1.1*H)
95+
paperParams(1,1)
96+
par(mar = c(3.7,3.7,0.8,0.8)); # beyond plot frame [b,l,t,r]
97+
par(mgp = c(2,0.8,0)); # placement of axis labels [1], numbers [2] and symbols[3]
98+
paperPlot(xlimit=c(2.8333,max(times2)), ylimit=c(0,100), xlabel='Time [hpi]', ylabel='% of Population', plotaxes=FALSE, xcol='black', ycol='black')
99+
resultsWT <- results1[1:3]
100+
lines(times2, lapply(resultsWT, fun, hashKey='XY % ++'), type='l', col='goldenrod3', lwd=2)
101+
lines(times2, lapply(resultsWT, fun, hashKey='XY % +-'), col = 'darkgreen', lwd = 2)
102+
lines(times2, lapply(resultsWT, fun, hashKey='XY % -+'), col = 'darkred', lwd = 2)
103+
lines(times2, lapply(resultsWT, fun, hashKey='XY % --'), col = 'black', lwd = 2)
104+
points(times2, lapply(resultsWT, fun, hashKey='XY % ++'), col='goldenrod3', lwd=2, cex=0.75)
105+
points(times2, lapply(resultsWT, fun, hashKey='XY % +-'), col = 'darkgreen', lwd = 2, cex=0.75)
106+
points(times2, lapply(resultsWT, fun, hashKey='XY % -+'), col = 'darkred', lwd = 2, cex=0.75)
107+
points(times2, lapply(resultsWT, fun, hashKey='XY % --'), col = 'black', lwd = 2, cex=0.75)
108+
axis(1, at=c(3,6,12,18))
109+
axis(2)
110+
box(col='black', lwd=1.5)
111+
legend('left',c('[+,+]','[+,-]','[-,+]','[-,-]'),lty=c(1,1,1,1),col=c('goldenrod3','darkgreen','darkred','black'), text.col=c('goldenrod3','darkgreen','darkred','black'), lwd=c(2,2,2), cex=0.85)
112+
dev.off()
113+
114+
######## Plot M51R ###############
115+
fun <- function(hashObj, hashKey)
116+
{
117+
return(hashObj[[hashKey]])
118+
}
119+
paperParams(1,1)
120+
pdf('M51R_FACSPercentages_Plot.pdf',width=1.1*W,height=1.1*H)
121+
paperParams(1,1)
122+
par(mar = c(3.7,3.7,0.8,0.8)); # beyond plot frame [b,l,t,r]
123+
par(mgp = c(2,0.8,0)); # placement of axis labels [1], numbers [2] and symbols[3]
124+
paperPlot(xlimit=c(2.8333,max(times2)), ylimit=c(0,100), xlabel='Time [hpi]', ylabel='% of Population', plotaxes=FALSE, xcol='black', ycol='black')
125+
resultsWT <- results1[4:6]
126+
lines(times2, lapply(resultsWT, fun, hashKey='XY % ++'), col='goldenrod3', lwd=2)
127+
lines(times2, lapply(resultsWT, fun, hashKey='XY % +-'), col = 'darkgreen', lwd = 2)
128+
lines(times2, lapply(resultsWT, fun, hashKey='XY % -+'), col = 'darkred', lwd = 2)
129+
lines(times2, lapply(resultsWT, fun, hashKey='XY % --'), col = 'black', lwd = 2)
130+
points(times2, lapply(resultsWT, fun, hashKey='XY % ++'), col='goldenrod3', lwd=2, cex=0.75)
131+
points(times2, lapply(resultsWT, fun, hashKey='XY % +-'), col = 'darkgreen', lwd = 2, cex=0.75)
132+
points(times2, lapply(resultsWT, fun, hashKey='XY % -+'), col = 'darkred', lwd = 2, cex=0.75)
133+
points(times2, lapply(resultsWT, fun, hashKey='XY % --'), col = 'black', lwd = 2, cex=0.75)
134+
axis(1, at=c(3,6,12,18))
135+
axis(2)
136+
box(col='black', lwd=1.5)
137+
legend('left',c('[+,+]','[+,-]','[-,+]','[-,-]'),lty=c(1,1,1,1),col=c('goldenrod3','darkgreen','darkred','black'), text.col=c('goldenrod3','darkgreen','darkred','black'), lwd=c(2,2,2), cex=0.85)
138+
dev.off()
139+
140+
######## Plot MOCK ###############
141+
fun <- function(hashObj, hashKey)
142+
{
143+
return(hashObj[[hashKey]])
144+
}
145+
paperParams(1,1)
146+
pdf('MOCK_FACSPercentages_Plot.pdf',width=1.1*W,height=1.1*H)
147+
paperParams(1,1)
148+
par(mar = c(3.7,3.7,0.8,0.8)); # beyond plot frame [b,l,t,r]
149+
par(mgp = c(2,0.8,0)); # placement of axis labels [1], numbers [2] and symbols[3]
150+
paperPlot(xlimit=c(2.8333,max(times2)), ylimit=c(0,100), xlabel='Time [hpi]', ylabel='% of Population', plotaxes=FALSE, xcol='black', ycol='black')
151+
resultsWT <- results1[7:9]
152+
lines(times2, lapply(resultsWT, fun, hashKey='XY % ++'), col='goldenrod3', lwd=2)
153+
lines(times2, lapply(resultsWT, fun, hashKey='XY % +-'), col = 'darkgreen', lwd = 2)
154+
lines(times2, lapply(resultsWT, fun, hashKey='XY % -+'), col = 'darkred', lwd = 2)
155+
lines(times2, lapply(resultsWT, fun, hashKey='XY % --'), col = 'black', lwd = 2)
156+
points(times2, lapply(resultsWT, fun, hashKey='XY % ++'), col='goldenrod3', lwd=2, cex=0.75)
157+
points(times2, lapply(resultsWT, fun, hashKey='XY % +-'), col = 'darkgreen', lwd = 2, cex=0.75)
158+
points(times2, lapply(resultsWT, fun, hashKey='XY % -+'), col = 'darkred', lwd = 2, cex=0.75)
159+
points(times2, lapply(resultsWT, fun, hashKey='XY % --'), col = 'black', lwd = 2, cex=0.75)
160+
axis(1, at=c(3,6,12,18))
161+
axis(2)
162+
box(col='black', lwd=1.5)
163+
legend('left',c('[+,+]','[+,-]','[-,+]','[-,-]'),lty=c(1,1,1,1),col=c('goldenrod3','darkgreen','darkred','black'), text.col=c('goldenrod3','darkgreen','darkred','black'), lwd=c(2,2,2), cex=0.85)
164+
dev.off()

0 commit comments

Comments
 (0)