dPHZ* PYO titration Dap vs. Dphys
09_11_18
knitr::opts_knit$set(root.dir = '/Users/scottsaunders/OneDrive - California Institute of Technology/Documents/Summer 2018/IDA_biofilms/09_11_18_DapEstimation')
knitr::opts_chunk$set(tidy.opts=list(width.cutoff=60),tidy=TRUE)
knitr::opts_chunk$set(echo = TRUE, fig.align="center")
suppressMessages(library('tidyverse'))
suppressMessages(library('cowplot'))
suppressMessages(library('broom'))
suppressMessages(library(modelr))
suppressMessages(library(viridis))
theme_set(theme_bw())Experiment
Now that I have repeated the dPHZ* titration, allowing the transferred biofilms to reach equilibrium, I wanted to try extracting actual diffusion coefficients. The simple idea is that if we can get a physical diffusion coefficient for PYO in the biofilm (Dphys) and we can get an apparent diffusion coefficient (Dap) for the electrons during the echem measurement then we can compare them. If physical diffusion of PYO is sufficient to explain the electrochemical signal then Dap ~ Dphys. However, if Dap >> Dphys that would imply there is another mechanism allowing electrons to diffuse.
Here I’m using already acquired data from 08/23 and 09/04.
Results
Estimate Dap from SWV vs. GC titration slope
In past notebooks I have shown that we can reproducibly get titration curves plotting GC current vs. SWV current. The slope of a linear fit of the data on that plot should yield Dap based on the following math:
C = concentration , D_ap = apparent diffusion coefficient , A = electrode area (SWV) , S = geometric factor (GC) , psi = SWV peak scalar , t_p = SWV pulse width , n= number of electrons per molecule , F = faraday’s constant, I=peak current
Slope of the linear fit of GC vs. SWV titration plot: \[slope = m = \frac{\Delta I_{gc}}{\Delta I_{swv}} \tag{eq. 1}\] Generator collector peak current: \[I_{gc} = nFSCD_{ap} \tag{eq. 2}\]
Square Wave peak current: \[I_{swv} = \frac{nFACD_{ap}^{1/2}}{\pi^{1/2}t_p^{1/2}} \psi \tag{eq. 3}\] Substituting into eq. 1: \[m = \frac{I_{gc}}{I_{swv}} = \frac{nFSCD_{ap}}{\frac{nFACD_{ap}^{1/2}}{\pi^{1/2}t_p^{1/2}} \psi} \tag{eq. 4}\] Simplifying: \[m = \frac{S\pi^{1/2}t_p^{1/2}D_{ap}^{1/2}}{A\psi} \tag{eq. 5}\]
\[m^2 = ({\frac{S\pi^{1/2}t_p^{1/2}D_{ap}^{1/2}}{A\psi}})^2 \tag{eq. 6}\]
\[m^2 = \frac{S^2\pi t_p D_{ap}}{A^2\psi^2} \tag{eq. 7}\] Rearranged for Dap: \[D_{ap} = \frac{m^2 A^2 \psi^2}{S^2\pi t_p} \tag{eq. 8}\]
Now we just need to process the data to get out the slope, and then we can calculate Dap given the other parameters.
Processed Titration Data from triplicate phz biofilms
# get already quantified GC and SWV data from experiments on
# 8/23 and 09/04/18 S1 = segment 1, S2 = segment 2 from the
# GC data
swvGCdataS1_Day1 <- read.csv("swvGCdataS1_08_23_18.csv")
swvGCdataS2_Day1 <- read.csv("swvGCdataS2_08_23_18.csv")
swvGCdataS1_Day2 <- read.csv("swvGCdataS1_09_04_18.csv")
swvGCdataS2_Day2 <- read.csv("swvGCdataS2_09_04_18.csv")
# Clean columns from data to prepare for merging both Days
swvGCdataS2_Day1 <- swvGCdataS2_Day1 %>% select(strain, reactor,
condition, PHZadded, PHZaddedInt, electrode, swvPYOsignal,
swvPYOsignalLast, gcPYOsignal) %>% mutate(swvElectrode = electrode) %>%
select(-electrode)
swvGCdataS2_Day2 <- swvGCdataS2_Day2 %>% select(strain, reactor,
condition, PHZadded, PHZaddedInt, electrode.x, swvPYOsignal,
gcPYOsignal) %>% mutate(swvElectrode = electrode.x, swvPYOsignalLast = NA,
reactor = "C") %>% mutate(condition = if_else(condition ==
"tran", "transfer", "soak")) %>% select(-electrode.x)
swvGCmultiDayS2 <- bind_rows(swvGCdataS2_Day1, swvGCdataS2_Day2)
ggplot(swvGCmultiDayS2 %>% filter(PHZaddedInt > 0 & swvElectrode ==
"i1"), aes(x = swvPYOsignal, y = gcPYOsignal, color = reactor,
fill = reactor)) + geom_line() + geom_point(shape = 21, color = "black") +
facet_wrap(~condition, scales = "free") + labs(title = "GC current vs. SWV current for Replicate dPHZ biofilms A,B, and C (S2 , i1)")
Seeing the data like this shows how consistent the soak data is. The data looks very linear and is very consistent across replicates. There is definitely more variability in the transfer steps, but it could be argued that the data is pretty linear…therefore it’s appropriate to perform a linear fit and examine the slopes.
Linear fits for triplicate titrations
Soak vs. Transfer
ggplot(swvGCmultiDayS2 %>% filter(PHZaddedInt > 0 & swvElectrode ==
"i1"), aes(x = swvPYOsignal, y = gcPYOsignal, linetype = condition,
color = reactor, fill = reactor)) + geom_smooth(method = "lm",
se = FALSE, aes(), fullrange = T) + geom_point(shape = 21,
color = "black") + labs(title = "Linear fits for soak step vs. transfer step - GC current vs. SWV current (S2 , i1)")
This plot shows that despite the variability in the transfer data, the slopes of the linear fits are always significantly less than the soak steps. This is reasonable and has been extremely consistent.
Transfer Only
ggplot(swvGCmultiDayS2 %>% filter(PHZaddedInt > 0 & swvElectrode ==
"i1" & condition == "transfer"), aes(x = swvPYOsignal, y = gcPYOsignal,
color = reactor, fill = reactor, linetype = condition)) +
geom_smooth(method = "lm", se = FALSE, fullrange = T) + geom_point(shape = 21,
color = "black") + facet_grid(~reactor) + labs(title = "Linear fits for only transfer step - GC current vs. SWV current (S2 , i1)")
Looking at the individual fits for the transfer titrations, it is clear that the fits are not perfect. In the future we may want to only fit parts of the dataset? The only reassuring thing is that the differences between the soak and transfer conditions swamp the variability between the transfer titrations…
# function to go with map() to fit lm model to each subset of
# data
multiDayBy <- swvGCmultiDayS2 %>% filter(swvElectrode == "i1" &
PHZaddedInt > 0) %>% group_by(condition, reactor) %>% nest()
lmModel <- function(df) {
lm(gcPYOsignal ~ swvPYOsignal, data = df, model = T)
}
# map the model function (above) to the nested data, then add
# predictions from each of those models to the appropriate
# dataset
multiDayBy <- multiDayBy %>% mutate(model = map(data, lmModel)) %>%
mutate(preds = map2(data, model, add_predictions)) %>% mutate(resids = map2(data,
model, add_residuals))
multiDayPredsS2 <- unnest(multiDayBy, preds, resids)
ggplot(multiDayPredsS2, aes(x = swvPYOsignal, y = gcPYOsignal,
color = condition)) + geom_point(aes(shape = reactor)) +
geom_line(aes(y = pred, linetype = reactor)) #+
coord_cartesian(ylim = c(0, 5e-08), xlim = c(0, 1.5e-07))
ggplot(multiDayPredsS2, aes(x = swvPYOsignal, y = gcPYOsignal,
color = condition)) + geom_point(aes(shape = reactor)) +
geom_line(aes(y = pred, linetype = reactor)) + facet_wrap(reactor ~
condition, scales = "free")
ggplot(multiDayPredsS2 %>% filter(condition == "transfer"), aes(x = swvPYOsignal,
y = gcPYOsignal, color = condition)) + geom_point(aes(shape = reactor)) +
geom_line(aes(y = pred, linetype = reactor)) + facet_wrap(~reactor,
scales = "free")lmsMultiDayS2 <- swvGCmultiDayS2 %>% filter(swvElectrode == "i1" &
PHZaddedInt > 0) %>% group_by(condition, reactor) %>% do(tidy(lm(gcPYOsignal ~
swvPYOsignal, .)))
# ggplot(lmsMultiDayS2 %>%
# filter(term=='swvPYOsignal'),aes(x=reactor,y=estimate,color=condition))+
# geom_point()+ylim(0,NA)
psi = 0.91
A = 0.039
S = 14.5
S = 20.4
tp = 0.0333
# pi
lmsMultiDayS2 <- lmsMultiDayS2 %>% filter(term == "swvPYOsignal") %>%
mutate(Dap = (estimate * A * psi)^2/(S^2 * (pi * tp))) %>%
group_by(condition) %>% mutate(mean = mean(Dap), Dcat = "Dap")
mean1 <- lmsMultiDayS2$mean[1]
mean2 <- lmsMultiDayS2$mean[4]
soakLabel <- paste("Dap = ", format(mean1, digits = 3), "cm^2/s")
tranLabel <- paste("Dap = ", format(mean2, digits = 3), "cm^2/s")
ggplot(lmsMultiDayS2, aes(x = condition, y = Dap, fill = reactor)) +
geom_col(position = "dodge", color = "black") + geom_point(aes(y = mean)) +
ylim(0, NA) + annotate("text", x = c(1.5, 2), y = c(mean1,
mean2 + 1.5e-06), label = c(soakLabel, tranLabel)) + labs(title = "Calculated Dap from titration slopes: Soak vs. Transfer")
This plot shows that the Dap measured from the titrations varies significantly between the soak and transfer condition. The Dap calculated for the soak step is around what we would expect for PYO physically diffusing in solution (~1 x 10^-5 cm^2/s).
\[D_{ap} = \frac{m^2 A^2 \psi^2}{S^2\pi t_p}\]
Parameters Used:
\(A= 0.039\) cm^2 (SWV Area)
\(S= 20.4\) cm (GC electrode geometry)
\(t_p= 0.0333\) sec (SWV pulse width)
\(\psi = 0.91\) (SWV peak scalar)
Note it should also be possible to propagate the uncertainty from the slope estimation to the calculated Dap…that way we could get some error bars on the above plot. In general, I think the uncertainty is pretty small though, because the slopes are all statistically significant.
Estimate Dphys for PYO, from SWV decay curves
Dilution of an initial mass with 1D diffusion next to a no flux boundary. \[C(t) = \frac{2M_0}{A_{xy}\sqrt{4\pi D_{phys}t}} \tag{eq. 1}\] Peak Square Wave current is given by: \[I_{swv} = \frac{nFA_{xy}C \sqrt{D_{ap}}}{\sqrt{\pi t_p}} \psi \tag{eq. 2}\] Concentration can be broken into mass / volume \[C = \frac{M}{V} = \frac{I_{swv}\sqrt{\pi t_p}}{nFA_{xy}\sqrt{D_{ap}}\psi} \tag{eq. 3}\] Therefore initial mass would be: \[M_0 = \frac{I_0\sqrt{\pi t_p}V}{nFA_{xy}\sqrt{D_{ap}}\psi} \tag{eq. 4}\] Substituting into eq. 1 gives: \[\frac{I_{swv}\sqrt{\pi t_p}}{nFA_{xy}\sqrt{D_{ap}}\psi} = \frac{2 \frac{I_0\sqrt{\pi t_p}V}{nFA_{xy}\sqrt{D_{ap}}\psi}}{A_{xy}\sqrt{4\pi D_{phys}t}} \tag{eq. 5}\] Simplified: \[I_{swv} = \frac{2I_0V}{A_{xy} \sqrt{4 \pi D_{phys}t}} \tag{eq. 6}\] The relevant volume should be the surface area of the electrode * the height of the biofilm? \[\frac{V}{A_{xy}}=z \tag{eq. 7}\] Below I plot this normalized decay curve: \[\frac{I_{swv}}{I_0} = \frac{z}{\sqrt{ \pi D_{phys}}}t^{-1/2} \tag{eq. 8}\] Fitting the normalized decay curve againt t^-(1/2) yields this slope: \[slope = m = \frac{z}{\sqrt{ \pi D_{phys}}} \tag{eq. 9}\] Rearranging gives Dphys: \[D_{phys} = \frac{z^2}{\pi m^2} \tag{eq. 10}\]
Note that this Dphys value is very dependent on the z value we pick. My intuition is that this z value should reflect the thickness of the biofilm…but I’m not totally sure. Even if it is, how accurately can we know that thickness…the biofilm is very heterogeneous.
files1 <- paste("swvDecay/reactorA", dir(path = "/Users/scottsaunders/OneDrive - California Institute of Technology/Documents/Summer 2018/IDA_biofilms/09_11_18_DapEstimation/swvDecay/reactorA",
pattern = "*.txt"), sep = "/")
df1 <- data_frame(filename = files1) %>% mutate(file_contents = map(filename,
~read_csv(., col_names = c("E", "i1", "N", "i2"))))
## Removed 25uM rep 8 from swv decay A because it had both
## negative and positive scans...
swvDecayA <- unnest(df1) %>% separate(filename, c("subDir", "subDir2",
"filename"), sep = "/") %>% separate(filename, c("PHZadded",
"echem", "rep"), sep = "_") %>% mutate(PHZaddedInt = as.integer(str_extract(PHZadded,
"^[0-9]+"))) %>% mutate(rep = as.integer(str_extract(rep,
"^[0-9]+"))) %>% mutate(reactor = "A") %>% select(PHZadded,
PHZaddedInt, reactor, rep, i1, i2, E) %>% gather(key = electrode,
value = current, i1, i2)
files2 <- paste("swvDecay/reactorB", dir(path = "/Users/scottsaunders/OneDrive - California Institute of Technology/Documents/Summer 2018/IDA_biofilms/09_11_18_DapEstimation/swvDecay/reactorB",
pattern = "*.txt"), sep = "/")
df2 <- data_frame(filename = files2) %>% mutate(file_contents = map(filename,
~read_csv(., col_names = c("E", "i1", "i2"))))
swvDecayB <- unnest(df2) %>% separate(filename, c("subDir", "subDir2",
"filename"), sep = "/") %>% separate(filename, c("PHZadded",
"condition", "echem", "rep"), sep = "_") %>% mutate(PHZaddedInt = as.integer(str_extract(PHZadded,
"^[0-9]+"))) %>% mutate(rep = as.integer(str_extract(rep,
"^[0-9]+"))) %>% mutate(reactor = "B") %>% select(PHZadded,
PHZaddedInt, reactor, rep, i1, i2, E) %>% gather(key = electrode,
value = current, i1, i2)
files3 <- paste("swvDecay/reactorC", dir(path = "/Users/scottsaunders/OneDrive - California Institute of Technology/Documents/Summer 2018/IDA_biofilms/09_11_18_DapEstimation/swvDecay/reactorC",
pattern = "*.txt"), sep = "/")
df3 <- data_frame(filename = files3) %>% mutate(file_contents = map(filename,
~read_csv(., col_names = c("E", "i1", "i2"))))
swvDecayC <- unnest(df3) %>% separate(filename, c("subDir", "subDir2",
"filename"), sep = "/") %>% separate(filename, c("strain",
"IDA", "reactor", "PHZadded", "molecule", "condition", "echem",
"rep"), sep = "_") %>% mutate(PHZaddedInt = as.integer(str_extract(PHZadded,
"^[0-9]+"))) %>% mutate(rep = as.integer(str_extract(rep,
"^[0-9]+"))) %>% mutate(reactor = "C") %>% select(PHZadded,
PHZaddedInt, reactor, rep, i1, i2, E) %>% gather(key = electrode,
value = current, i1, i2)
# Combine SWV decays from replicates A, B and C
swvDecayMulti <- bind_rows(swvDecayA, swvDecayB, swvDecayC) %>%
arrange(reactor, PHZaddedInt, E)
# ggplot(swvDecayMulti %>%
# filter(electrode=='i1'),aes(x=E,y=current,color=rep))+geom_path(aes(group=rep))+
# scale_color_viridis()+scale_x_reverse()+
# facet_grid(reactor~PHZaddedInt,scales='free')+
# labs(title='SWV signal decays in transfer reactor (after
# each soak) for biofilms A,B & C')swvPYOmaxMulti <- swvDecayMulti %>% filter(E >= -0.3 & E <= -0.2) %>%
group_by(reactor, PHZadded, PHZaddedInt, rep, electrode) %>%
mutate(PYOmax = max(current)) %>% filter(current == PYOmax)
swvPYOminMulti <- swvDecayMulti %>% filter(E >= -0.4 & E <= -0.1) %>%
group_by(reactor, PHZadded, PHZaddedInt, rep, electrode) %>%
mutate(PYOmin = min(current)) %>% filter(current == PYOmin)
ggplot(swvDecayMulti %>% filter(electrode == "i1"), aes(x = E,
y = current, color = rep)) + geom_path(aes(group = rep)) +
geom_point(data = swvPYOmaxMulti %>% filter(electrode ==
"i1"), aes(x = E, y = PYOmax), color = "red", size = 0.1) +
geom_point(data = swvPYOminMulti %>% filter(electrode ==
"i1"), aes(x = E, y = PYOmin), color = "cyan", size = 0.1) +
facet_grid(reactor ~ PHZaddedInt, scales = "free") + scale_color_viridis() +
labs(title = "raw SWV signals decay in transfer reactor (after each soak) for biofilms A,B & C")
# ggplot(swvDecayMulti %>%
# filter(electrode=='i1'),aes(x=E,y=current,color=rep))+geom_path(aes(group=rep))+
# geom_point(data=swvPYOmaxMulti %>%
# filter(electrode=='i1'),aes(x=E,y=PYOmax),color='red',size=0.1)+
# geom_point(data=swvPYOminMulti %>%
# filter(electrode=='i1'),aes(x=E,y=PYOmin),color='cyan',size=0.1)+
# facet_wrap(reactor~PHZaddedInt,scales='free')+labs(title='')+
# scale_color_viridis()We can find the min max points in all of the SWVs and they look reasonably consistent - signal decreases with decreasing soak concentration.
swvPYOsigMulti <- left_join(swvPYOmaxMulti, swvPYOminMulti, by = c("reactor",
"electrode", "PHZadded", "PHZaddedInt", "rep")) %>% ungroup() %>%
mutate(swvPYOsignal = PYOmax - PYOmin) %>% arrange(reactor,
electrode, rep) %>% mutate(t = (rep - 0.5) * 105, t0.5 = t^-0.5,
t1 = t^-1)
ggplot(swvPYOsigMulti %>% filter(electrode == "i1"), aes(x = t,
y = swvPYOsignal, color = PHZaddedInt)) + geom_point() +
geom_line(aes(group = PHZaddedInt)) + facet_grid(~reactor,
scales = "free") + scale_color_viridis() + labs(title = "Quantified SWV decay curves")
Decay curves look pretty good. It is interesting that replicate C always starts at significantly higher concentrations/signals…not sure why that is. One slightly concerning thing is that for replicate C, transfers 150 and 200uM, I might not have waited long enough for them to reach equilibrium. This might be the reason that they show abnormally high GC/SWV values in the GC vs. SWV section above.
Decay curves are linear vs. t^-(1/2)
Curves
ggplot(swvPYOsigMulti %>% filter(electrode == "i1"), aes(x = t^-0.5,
y = swvPYOsignal, color = PHZaddedInt)) + geom_point() +
geom_line(aes(group = PHZaddedInt)) + facet_wrap(~reactor,
scales = "free") + scale_color_viridis() + labs(title = "SWV decay curves vs. t^-(1/2)")
Decay curves look pretty linear, althouh the highest points on each plot seem to throw off the trend slightly.
Linear Fits
ggplot(swvPYOsigMulti %>% filter(electrode == "i1" & t > 105),
aes(x = t^-0.5, y = swvPYOsignal, color = PHZaddedInt)) +
geom_point() + geom_smooth(aes(group = PHZaddedInt), method = "lm") +
facet_wrap(~reactor, scales = "free") + scale_color_viridis() +
labs(title = "SWV decay Linear Fit vs. t^-(1/2) (excluding first point)")
Exluding the first titration points, the decay curves look extremely linear - Fits look very good.
Normalized Decay Curves
Normalizing the Decay curves by the initial signal and plotting against t^-(1/2) allows us to compare the slopes that we will use to calculate Dphys. Below I am showing one example from the 50uM transfer condition.
Normalized vs. t
maxMulti <- swvPYOsigMulti %>% filter(t > 105) %>% group_by(reactor,
PHZadded, PHZaddedInt, electrode) %>% mutate(I0 = max(swvPYOsignal))
ggplot(maxMulti %>% filter(electrode == "i1" & PHZaddedInt ==
50), aes(x = t, y = swvPYOsignal/I0, color = reactor, fill = reactor)) +
geom_line(aes(group = reactor)) + geom_point(shape = 21,
color = "black") + labs(title = "Normalized decay curves (50uM soak)")
Linear fit of normalized vs. t^-(1/2)
ggplot(maxMulti %>% filter(electrode == "i1" & PHZaddedInt ==
50), aes(x = t^-0.5, y = swvPYOsignal/I0, color = reactor,
fill = reactor)) + geom_smooth(aes(group = reactor), method = "lm",
se = F) + geom_point(shape = 21, color = "black") + labs(title = "Linear fit for normalized decay curves vs. t^-(1/2) (50uM soak)")
Plotting against this time scale shows that the linear fits look reasonable and that the slopes look pretty similar between replicates.
Fitting the whole dataset in the same fashion and extracting the slope values allows us to estimate Dphys using this equation: \[D_{phys} = \frac{z^2}{\pi m^2}\]
z = 10^-3 #10um
lmMulti <- maxMulti %>% group_by(reactor, PHZadded, PHZaddedInt,
electrode) %>% do(tidy(lm(swvPYOsignal/I0 ~ t0.5, .))) %>%
mutate(Dphys = (z^2)/(pi * estimate^2), Dcat = "Dphys") %>%
mutate(condition = "transfer")
z100 = 10^-2 #100um
lmMultiZ100 <- maxMulti %>% group_by(reactor, PHZadded, PHZaddedInt,
electrode) %>% do(tidy(lm(swvPYOsignal/I0 ~ t0.5, .))) %>%
mutate(Dphys = (z100^2)/(pi * estimate^2), Dcat = "Dphys") %>%
mutate(condition = "transfer")
ggplot(lmMulti %>% filter(term == "t0.5" & PHZaddedInt > 5),
aes(x = factor(PHZaddedInt), y = Dphys, fill = reactor)) +
geom_bar(stat = "summary", fun.y = "mean", position = "dodge",
color = "black") + geom_point(aes(shape = reactor), position = position_dodge(width = 0.75)) +
labs(y = "Calculated Dphys (cm^2)/s", x = "PYO added during soak (uM)",
title = "Calculated physical diffusion coefficients")
For \(z=0.001\) cm (\(z=10\) um)
It is reassuring that we get values in the same order of magnitude for all the different concentrations and replicates and electrode (each SWV i1 and i2 are shown as the black dots). Each replicate shows slightly different trends (e.g. C shows lower D at higher [PYO]), but I think the take away should be that it’s pretty consistent.
Surprisingly, these seem like somewhat reasonable values. The estimated diffusion time for a molecule with D=2x10^-9 cm^2/s to move 25um is about 26 min. We wait about an hour to make sure the PYO is at equilibrium, but that sounds like the right order of magnitude…
Meanwhile, for the measured Dap during soak (D~8x10^-6 cm^2/s) the diffusion time across 25um would be 0.39 seconds. That also seems reasonable given that we almost instantly see no PCA in the transferred biofilm.
http://www.physiologyweb.com/calculators/diffusion_time_calculator.html
Comparing calculated Dap and Dphys
Let’s compare on the same plot our Dap and Dphys estimates. Note these plots are log scale. Remember that the Dphys estimate is highly dependent on the value for z, so let’s look at a z of 10um and 100um.
z=10um
lmMultiSum <- lmMulti %>% filter(term == "t0.5") %>% group_by(reactor) %>%
summarise(Dap = mean(Dphys)) %>% mutate(Dcat = "Dphys", condition = "transfer Dphys")
dCalcCompare <- bind_rows(lmMultiSum, lmsMultiDayS2 %>% select(reactor,
Dap) %>% ungroup() %>% mutate(Dcat = "Dap", condition = fct_recode(condition,
soakDap = "soak", transferDap = "transfer"))) %>% mutate(condition = fct_relevel(condition,
"soakDap", "transferDap", "transfer Dphys"))
ggplot(dCalcCompare, aes(x = condition, y = Dap, color = reactor)) +
geom_boxplot(aes(group = condition)) + geom_jitter(height = 0,
width = 0.25, size = 3) + scale_y_log10() + labs(x = "",
y = "Calculated D (cm^2/s)", title = "Comparing Calculated Diffusion Coefficients \nDap (soak and transfer) and Dphys | z=10um")
z=100um
lmMultiZ100Sum <- lmMultiZ100 %>% filter(term == "t0.5") %>%
group_by(reactor) %>% summarise(Dap = mean(Dphys)) %>% mutate(Dcat = "Dphys",
condition = "transfer Dphys")
dCalcCompareZ100 <- bind_rows(lmMultiZ100Sum, lmsMultiDayS2 %>%
select(reactor, Dap) %>% ungroup() %>% mutate(Dcat = "Dap",
condition = fct_recode(condition, soakDap = "soak", transferDap = "transfer"))) %>%
mutate(condition = fct_relevel(condition, "soakDap", "transferDap",
"transfer Dphys"))
ggplot(dCalcCompareZ100, aes(x = condition, y = Dap, color = reactor)) +
geom_boxplot(aes(group = condition)) + geom_jitter(height = 0,
width = 0.25, size = 3) + scale_y_log10() + labs(x = "",
y = "Calculated D (cm^2/s)", title = "Comparing Calculated Diffusion Coefficients \nDap (soak and transfer) and Dphys | z=100um")
Conclusions
- Dap and Dphys estimates yield reasonable values that show expected and unexpected differences. We expected Dap to be higher for the soak step than the transfer step, but the difference between Dphys and Dap was unexpected (by me…).
- In order to believe these results we need to:
- believe our values for z/V/A
- Show that the Dphys method is capable of recovering realistic diffusion coefficients in vitro.
- Test if the Dap values derived from potential step match this GC vs. SWV method.