Overview

Here we investigate temporal patterns in floral volatile emissions in Schiedea hookeri, Schiedea kaalae, and their hybrid. Prior GC-MS analysis and the human nose suggest differences in total volatiles output and composition between broad day and evening sampling periods. PTR-MS allows quantification of these volatile emission rates in real time. This report includes:

  • photosynthetic parameters inferred from carbon dioxide and water fluxes using an infrared gas analyzer (IRGA)
  • timeseries of individual ion peaks: small molecules, known top compounds from each species, other ions of interest
  • a clustering analysis of timeseries for all peaks to group ions by their emission patterns

The PTR-MS produces ions via proton transfer from hydronium ions (H3O+). The ions are transfered from the source into a drift tube and a time-of-flight mass spectrometer. Analytes undergo soft chemical ionization and do not intiailly fragment, but due to the electric field, ions may fragment in the drift tube. Mono- and sesquiterpenes are known to fragment depending on the strength of this field (Faiola et al., in prep).

A PTR-Quadrupole (actual instrument used is a PTR-TOF)

A PTR-Quadrupole (actual instrument used is a PTR-TOF)

First moth visits

vst <- read.table("firstvisits.csv", header=T, sep="\t", colClasses=c(firstmoth="POSIXct"))

library(maptools)
ekahanui <- matrix(c(-158.088214, 21.440004), nrow=1)
vst$sunsets <- do.call("c", lapply(lapply(vst$firstmoth, sunriset, crds=ekahanui, direction="sunset", POSIXct.out=T), function(x) x$time))
attr(vst$sunsets, "tz") <- Sys.timezone()
vst$sunsets  <- as.POSIXct(format(vst$sunsets, tz="US/Hawaii", usetz=TRUE))
vst$firstmothsunset <-vst$firstmoth - vst$sunsets

vst$solarnoons <- do.call("c", lapply(lapply(vst$firstmoth, solarnoon, crds=ekahanui, POSIXct.out=T), function(x) x$time))
attr(vst$solarnoons, "tz") <- Sys.timezone()
vst$solarnoons  <- as.POSIXct(format(vst$solarnoons, tz="US/Hawaii", usetz=TRUE))
vst$firstmothnoon <- vst$firstmoth - vst$solarnoons

kable(vst)
firstmoth sunsets firstmothsunset solarnoons firstmothnoon
2014-03-03 17:48:00 2014-03-03 18:37:31 -49.51667 mins 2014-03-03 12:44:11 5.063611 hours
2014-03-04 18:10:00 2014-03-04 18:37:54 -27.90000 mins 2014-03-04 12:43:58 5.433889 hours
2014-03-05 18:00:00 2014-03-05 18:38:17 -38.28333 mins 2014-03-05 12:43:44 5.271111 hours
2014-07-03 17:58:00 2014-07-03 19:18:48 -80.80000 mins 2014-07-03 12:36:39 5.355833 hours
2015-02-16 17:51:00 2015-02-16 18:30:45 -39.75000 mins 2015-02-16 12:46:25 5.076389 hours
2015-02-17 17:50:00 2015-02-17 18:31:15 -41.25000 mins 2015-02-17 12:46:21 5.060833 hours
2015-02-18 17:56:00 2015-02-18 18:31:44 -35.73333 mins 2015-02-18 12:46:16 5.162222 hours
2015-02-19 17:48:00 2015-02-19 18:32:13 -44.21667 mins 2015-02-19 12:46:10 5.030556 hours
2015-02-21 17:52:00 2015-02-21 18:33:09 -41.15000 mins 2015-02-21 12:45:57 5.100833 hours
2016-02-09 17:47:00 2016-02-09 18:26:54 -39.90000 mins 2016-02-09 12:46:34 5.007222 hours
qplot(as.numeric(vst$firstmothsunset), xlim=c(-100,0), binwidth=10, xlab="Minutes until sunset",ylab="Count (n=10)", main="First Pseudoschankia brevipalpis visits to S. kaalae flowers at ‘Ēkahanui Gulch") +theme_bw()

Experimental

sp <- read.delim("Schiedea_Fluxtron_20171106_sync.csv")
rownames(sp) <- paste0(sp$Species, sp$Experiment, sp$Cuvette)
nexp <- length(unique(sp$Experiment))
exper <- !duplicated(sp$Experiment)
spcol <- setNames(c("grey20",brewer.pal(8, "Set1"), brewer.pal(3, "Dark2"))[c(1,3,2,6,5,8,9,4,10)], rownames(sp))
cscale <- scale_color_manual("Species",values=spcol)
cguide <- guides(color = guide_legend(override.aes = list(size=3)))

tzone <- "Etc/GMT+8"
tfrmt <- "%Y-%m-%d %H:%M:%S"
sp$start <- as.POSIXct(strptime(sp$start, format=tfrmt, tz=tzone))
sp$lightoff <- as.POSIXct(strptime(sp$lightoff, format=tfrmt, tz=tzone))
tp <- mapply(seq.POSIXt, from=sp$start, length.out=sp$steps, MoreArgs=list(by="sec")) # changed starts in file to be the same date!
lights <- mapply(seq.POSIXt, from=sp$lightoff[exper], length.out=sp$changes[exper], MoreArgs=list(by=12 * 60 * 60))
tscale <- lapply(lights, scale_x_datetime, name="Time", date_minor_breaks="1 hour", date_labels="%H:%M",expand = c(0, 0))

rects <- lapply(1:nexp, function(i) {data.frame(xstart = c(sp$start[exper][i], lights[[i]]), xend = c(lights[[i]], tail(tp[exper][[i]], 1)))})
rects <- lapply(rects, function(x) cbind(x,Period=rep_len(c("day","night"), nrow(x)),stringsAsFactors=F))
grect <- lapply(rects, function(x) geom_rect(data=x, mapping=aes(xmin = xstart, xmax = xend, ymin = -Inf, ymax = Inf, fill = Period), alpha = 0.4, inherit.aes = FALSE))

greenhouse <- matrix(c(-117.8475746, 33.6472914), nrow=1)
sp$sunrise <- sapply(tp, function(x) {sunriset(crds=greenhouse, dateTime=x[1]-24*60*60, direction="sunrise", POSIXct.out=T)$time} )
sp$sunset <- sapply(tp, function(x) {sunriset(crds=greenhouse, dateTime=x[1]-24*60*60, direction="sunset", POSIXct.out=T)$time} )
attributes(sp$sunrise) <- attributes(sp$sunset) <- attributes(tp[[1]])

Three plants with open flowers were selected from greenhouse (one of each type). Plants were all watered prior to sampling and placed in a photochamber set to 23 C and with 60% relative humidity. During the photoperiod, plants received about 450 µmol/(m2 s) PAR (see table) at the top of the inflorescence from LED lights (c.f. direct sunlight = 2,000 µmol/(m2 s) PAR).

Flowering inflorescences from each plant were supported with a stake and placed in “cuvettes” made of nylon-6 oven bags (dimensions 40.6 cm x 44.4 cm). An empty bag was used as a blank (cuvette #2). Each cuvette recieved dry zero air at the bottom of the bag from a zero-air generator (critical orifice Reynolds number 33, 25 psi, 12 L/min of air combined to the 4 cuvettes). Air was sampled from the top of the bag. Each cuvette contailed a thermocouple to measure temperature in the bag. The bag was sealed with a twist-tie around the pedicel, input and sampling lines, and thermocouple.

All times are PST (UTC-8). The experiment spanned hours, from mbient sunrise and sunset times the day prior were

The Fluxtron switched gas flow to the PTR-MS between each cuvette every 4 min, and flow to the IRGA between the blank cuvette and sampled cuvette every 20 s. The IRGA measured CO2 concentration and temperature. The Ionicon PTR-TOF held the settings below for volatiles analysis.

Plants

kable(sp[,1:14], row.names=F)
Plant identity, and infloresence compostion. Flowers counted after experiment ended, at 2017-08-11 09:40 for Exp 0. Plants watered at 2017-08-09 ~10:30 and 2017-08-11 08:10 fro Exp 0.
Batch Experiment Cuvette Species Population ID Sex Leaves Inflor Mphase Fphase Closed Buds Flrs
2017-08-09 0 1 KH 3587-7 x 879-10-1 189-3 H 2 1 6 29 29 50 35
2017-08-09 0 3 hookeri 891 3-11 H 8 8 4 24 153 30 28
2017-08-09 0 4 kaalae 3587-A x 3587-14 229-3 H 0 3 15 39 114 83 54
2017-11-06 1 1 kaalae 904-2 x 904-5 158-3 H 0 1 2 3 203 8 5
2017-11-06 1 3 kaalae 904-2 x 904-5 158-5 H 0 1 1 51 136 89 52
2017-11-06 1 4 kaalae 904-5 x 904-3 184-18 H 0 2 17 46 377 220 63
2017-11-06 2 1 kaalae 3587 6-1 H 0 1 0 18 290 8 18
2017-11-06 2 3 hookeri 866 2 H 12 4 3 16 51 17 19
2017-11-06 2 4 hookeri 879 3-8 H 0 1 0 9 32 11 9

PTR-MS parameters

PTR-MS parameters. TPS settings of 2017-04-07. *Inlet flow set to 150 scmm until 2017-08-09 13:10.
Parameter Value
MCP 2100 V
inlet temperature 60 C
inlet flow 140 sccm*
H2O flow 5 sccm
ion source current 3 mA
drift tube temperature 60 C
drift tube pressure 2.8 mbar
Us 150 V
Uso 80 V
Udrift 600 V
Ufunnel 30 V
Uout 7 V
Urf 16 V

Photosynthesis parameters

Photosynthesis parameters calculated with the vonCaemmerer and Farquhar (1981) model. Results of fx FT_PS_LiCor_equations() on each of the rows of the FT waves.
Variable Units Description
Anet (mol CO2)/(m2 s) net CO2 assimilation rate
E (mol H2O)/(m2 s) water transpiration rate
Gsw (mol H2O)/(m2 s) stomatal conductance to water vapor
Gbw (mol H2O)//(m2 s) boundary layer conductance to water vapor
Gtw (mol H2O)/(m2 s) total conductance to water vapor
Ci (mol CO2)/(mol air) intercellular CO2 concentration
Gtc (mol CO2)/(m2 s) total conductance to CO2
Ue mol/s standardized flow of air entering the cuvette
SVP kPa saturated water vapor pressure at Tleaf
Wl (mol H2O)/(mol air) molar concentration of water vapor within the leaf

Environmental

Other data recorded by the fluxtron: * Flow to each cuvette (SLPM) * PAR light intensity (umol/m2/s) * CO2 concentration (ppm) * temperature in each cuvette (K)

flow.files = setNames(lapply(levels(sp$name), function(x) { all <- paste0(x, "/", dir(paste0(path, x))); paste0(path,all[grepl("Fluxtron2",all)]) }), levels(sp$name))
flow.cols <- colnames(fread(flow.files[[1]][1], header=T, sep=",", nrows=0))
read.flow <- function(file) {
fread(file, header=F, sep=",", skip=2, colClasses="numeric", col.names=flow.cols)
}
flow <- lapply(flow.files, function(x) { rbindlist(lapply(x,read.flow)) })
makeblocks <- function(x) { cbind(x, block=rep(1:nrow(x), each = 60, len=nrow(x))) }
flow <- lapply(flow, makeblocks)
flow <- lapply(flow, avgblocks)
conv.tzone <- function(x) { x$tp <- as.POSIXct(x$TS_Labview_UTC,  origin = "1904-01-01", tz = "GMT")+1; x } 
flow <- lapply(flow, conv.tzone)
flow <- lapply(flow, fixtz)
aligndays <- lapply(flow, function(x) { floor(difftime(x$tp[1], flow[[1]]$tp[1], "days")) })#assumes that the first file is the start day of the experiment
aligndates <- function(x, i) { x$tp <- x$tp - as.numeric(aligndays[i])*24*60*60; x}
flow <- mapply(aligndates, x=flow, i=1:3, SIMPLIFY=FALSE)
flow.plot <- function(fc) {
flow.plotvars <-   setNames(flow.cols[fc], flow.cols[fc])
mflow <- lapply(flow, melt, id.vars="tp")
mmflow <- do.call(rbind, lapply(names(mflow), function(x){data.frame(id=x, mflow[[x]])}))
mmflow <- mmflow[mmflow$variable %in% names(flow.plotvars),]
ggplot(mmflow, aes(x=tp, y=value, color=id)) + facet_wrap(~ variable, scales="free_y", ncol=1, labeller=labeller(variable=flow.plotvars)) + cguide  + scale_fill_manual(values=c("day"="grey60", "night"="grey40")) + geom_line(size=0.5, alpha=1) + tscale[[3]] + grect[[3]]
}

flow.plot(32:39) + ylim(1.4,2.8)#flows

flow.plot(28:31) + ylim(293,300) #temps

flow.plot(40:41)

Results

IRGA

Photosynthesis plots

irga.files = setNames(paste0(path, sp$name, "/FTirgaFluxc", sp$Cuvette, ".txt"), rownames(sp))
read.irga <- function(file, timepoints) {
  fixtz(medianblocks(blockjumps(na.omit(cbind(fread(file, header=T, sep="\t", skip=2, colClasses="numeric"), tp=timepoints)))))
}
irga <- mapply(read.irga, file=irga.files, timepoints=tp, SIMPLIFY=F, USE.NAMES=T)

#normailze by number of flowers
irga.normvars <- c("Anet","E","Gsw","Gbw","Gtw","Gtc")
irga.normalize <- function(x,flrs){x[,irga.normvars] <- x[,irga.normvars]/flrs; return(x)}
irgan <- mapply(irga.normalize, irga, sp$Flrs, SIMPLIFY=F)

irga.plotvars <-   c("Anet"="Anet: net CO2 assimilation rate",
                 "E"="E: water transpiration rate",
                 "Ue"="Ue: standardized flow of air entering the cuvette")
mirga <- lapply(irgan, melt, id.vars="tp")
mmirga <- do.call(rbind, lapply(names(mirga), function(x){data.frame(id=x, mirga[[x]])}))
mmirga <- mmirga[mmirga$variable %in% names(irga.plotvars),]

Calculated photosynthetic parameters Anet and E normalized by the number of flowers, and the flow of air entering the cuvette Ue.

irgaplot <- ggplot(mmirga, aes(x=tp, y=value, color=id)) + facet_wrap(~ variable, scales="free_y", ncol=1, labeller=labeller(variable=irga.plotvars)) + cscale + cguide + scale_fill_manual(values=c("day"="grey60", "night"="grey40"))
irgaplot + tscale[[3]] + grect[[3]] + geom_line(size=0.5, alpha=1)

#for(i in 1:nexp) {
#ip <- irgaplot %+% subset(mmirga, id %in% rownames(sp[sp$Experiment==i-1,])) + tscale[[i]] + grect[[i]] + geom_line(size=0.5, alpha=1)
#plot(ip)
#}

Photosynthesis discussion

Exp 0. Increase in noise at 4 am PST (Why 3 am in figure?) 17/08/10 . The CO2 concentration started to oscillate, with a markedly serrated pattern (compare the left-hand side of the graph to the right-hand side (black line is the signal from all cuvettes, the blue line is the blank coming from c2, and the red is the signal from c4).

Anet : net CO2 assimilation rate

  • Infloresences were net photosynthetic during the day, emitted CO<sub2~ through respiration at night
  • After normalizing by number of open flowers, Anet is similar
  • Responses <30 min

E: water transpiration rate

  • Water transpiration was higher during the day than night - opened stomata?
  • Responses <30 min

Ue: standardized flow of air entering the cuvette

  • S. kaalae flow is 30% higher than other two cuvettes - is this accounted for?

PTR-MS ion timecourses

ncps.files = setNames(paste0(path, sp$name, "/FTncpsFluxc", sp$Cuvette, ".txt"), rownames(sp))
read.ncps <- function(file, timepoints) {
  fixtz(medianblocks(blockjumps(na.omit(cbind(Filter(function(x) !all(is.na(x)), fread(file, header=T, sep="\t", skip=2, colClasses="numeric")), tp=timepoints)))))
}
#ncps <- mapply(read.ncps, file=ncps.files, timepoints=tp, SIMPLIFY=F, USE.NAMES=T)
load("ncps.Rdata")

library(lubridate)
remove.pib <- function(x) { 
  dec <- hour(x$tp)+minute(x$tp)/60
  x$pib <- FALSE#(dec > 16.9 & dec < 17.25) | (dec > 4.9 & dec < 5.2)
  x
}
ncps <- lapply(ncps, remove.pib)

#Weird blips
#ncps[["kaalae04"]]$pib[ncps[["kaalae04"]]$tp == as.POSIXct("2017-08-09 22:20:32.5", tz=tzone)] <- TRUE
#ncps[["kaalae11"]]$pib[ncps[["kaalae11"]]$tp == as.POSIXct("2017-08-12 04:47:46.5", tz=tzone)] <- TRUE
#ncps[["kaalae13"]]$pib[ncps[["kaalae13"]]$tp == as.POSIXct("2017-08-11 06:31:49.5", tz=tzone)] <- TRUE
#ncps[["kaalae14"]]$pib[ncps[["kaalae14"]]$tp == as.POSIXct("2017-08-09 15:55:53.5", tz=tzone)] <- TRUE
#ncps[["kaalae14"]]$pib[ncps[["kaalae14"]]$tp == as.POSIXct("2017-08-12 04:43:46.5", tz=tzone)] <- TRUE

noncount <- c("tp","block", "Group.1", "pib")
ncps.normalize <- function(x,flrs){x[,!colnames(x) %in% noncount] <- x[,!colnames(x) %in% noncount]/flrs; return(x)}
nncps <- mapply(ncps.normalize, ncps, sp$Closed+sp$Flrs, SIMPLIFY=F) #Norm by closed and open flowers - makes a difference for scales and thresholds

mncps <- lapply(nncps, melt, id.vars=c("tp","pib"))
mncps <- do.call(rbind, lapply(rownames(sp), function(x){data.frame(id=x, mncps[[x]])}))
#mncps$id <- factor(mncps$id,levels(mncps$id)[c(2,3,1)])
mdncps <- mncps
mdncps$day <- format(mncps$tp,"%d")
mdncps$tp <- as.POSIXct(paste("2017-08-10 ", format(mncps$tp,"%H:%M:%S")), tz=tzone)
frmla <- read.table("formulas.csv", sep="\t", header=T, stringsAsFactors=F, colClasses=c("factor"="factor", "cycle"="factor"))

Top volatiles

lineplot <- function(plotpeaks) {
ggplot(mncps[mncps$variable %in% names(plotpeaks) & !mncps$pib,], mapping=aes(x=tp, y=value, color=id)) + facet_wrap(~ variable, scales="free_y", ncol=1, labeller=labeller(variable=plotpeaks)) + tscale[[3]] + cscale + cguide + scale_fill_manual(values=c("day"="grey90", "night"="grey80")) + grect[[3]] + geom_line(size=0.5, alpha=1) }

dailyplot <- function(plotpeaks) {
ggplot(mdncps[mdncps$variable %in% names(plotpeaks) & !mncps$pib,], mapping=aes(x=tp, y=value, color=id, linetype=day)) + scale_linetype_manual(name="Day", values = c("longdash", "solid", "twodash","dotted", "dotdash")) + facet_wrap(~ variable, scales="free_y", ncol=1, labeller=labeller(variable=plotpeaks)) + tscale[[1]] + cscale + cguide + scale_fill_manual(values=c("day"="grey60", "night"="grey40")) + grect[[1]] + geom_line(size=0.5, alpha=1) }

HOOKpeaks <- c("129.127"="129: 1-octen-3-ol / 3-octanone / 2,3-heptanedione : HOOK 29.3%, KAAL 13.8%, KAHO 34.7%",
               "101.096"="101: 3-hexen-1-ol : HOOK 27.9%, KAAL 0.2%, KAHO 6.3%",
               "163.123"="163: nearC11H14O unknown : HOOK 15.3%, KAAL 0.0%, KAHO 0.7%",
               "107.049"="107: benzaldehyde : HOOK 5.2%, KAAL 0.9%, KAHO 0.7%") 
KAALpeaks <- c("169.122"="169: linalool ketone (pyranoid): HOOK 1.3%, KAAL 39.3%, KAHO 37.8%",
               "171.138"="171: linalool oxide (pyranoid): HOOK 0.2%, KAAL 17.4%, KAHO 3.6%",
               "243.159"="243: Linalool oxide (furanoid) ethyl carbonate: HOOK 0.2%, KAAL 13.9%, KAHO 3.0%",
               "121.065"="121: phenylacetaldehyde : HOOK 0.4%, KAAL 5.9%, KAHO 0.2%",
               "137.132"="137: α-phellandrene / monoterpenes : HOOK 0.1%, KAAL 2.6%, KAHO 0.4%")

MYSTpeaks <- c("163.123"="163: Unknown benzenoid","92.062"="92: Unknown benzenoid fragment 1","120.081"="120 Unknown benzenoid fragment 2")

lineplot(HOOKpeaks)

lineplot(KAALpeaks)

lineplot(MYSTpeaks)

Stacked plots

dailyplot(HOOKpeaks)

dailyplot(KAALpeaks)

All ions

ncpsp <- lapply(ncps, function(x) {x[,!names(x) %in% noncount]}) 
frmlas <- setNames(trimws(substr(paste(frmla$mass, frmla$factor, frmla$cycle, frmla$compound, frmla$common), 1, 40)), colnames(ncpsp[[1]]))
rownames(frmla) <- colnames(ncpsp[[1]])

ncps_plot <- function(pl) {
trange <- !ncps[[pl]]$pib
prange <- 31:ncol(ncpsp[[pl]])
ncpscut <-  ncpsp[[pl]][trange,prange]
ncpscut <- ncpscut[,colMax(ncpscut)>threshold]
ncpscut[ncpscut<0] <- 0
ncpscut <- log(ncpscut+1)
ncpscut.s <- divmax(ncpscut)
ncpscut.st <- t(as.matrix(ncpscut.s)) 
ncpscut.st3 <- ncpscut.st

ncpscut.tp <- ncps[[pl]]$tp[trange]

ncpscut.frmlas <- frmlas[colnames(ncpscut)]
ncpscut.frmlas <- aggregate(ncpscut.frmlas~round(mass), paste, collapse="/", data=frmla[colnames(ncpscut),])

wmass <- factor(round(frmla[colnames(ncpscut),"mass"]))
ncpscut.st <- as.data.frame(t(ncpscut.s))
ncpscut.st <- aggregate(. ~ wmass, mean, data=ncpscut.st)
rownames(ncpscut.st) <- ncpscut.st$wmass
ncpscut.st$wmass <- NULL

ncpscut.t <- as.data.frame(t(ncpscut))
ncpscut.t <- aggregate(. ~ wmass, mean, data=ncpscut.t)
rownames(ncpscut.t) <- ncpscut.t$wmass
ncpscut.t$wmass <- NULL

colside_max <- range01(rowMax(ncpscut.t)^(1/16))
sidecols <- colorRampPalette(c("white", "purple"))(512)[round(colside_max*511)+1]
xhclust <- function(data) { hclust(data, method="mcquitty") }
xdist <- function(data) { dist(data, method="DTW") }

#heatmap(as.matrix(ncpscut.st), Colv=NA, col=topo.colors(512), hclustfun=xhclust, distfun=xdist, labRow=rownames(ncpscut.st), labCol=NA, RowSideColors=sidecols, margins=c(5,18))
  
#rownames(ncpscut.st)<-ncpscut.frmlas[,2]
colnames(ncpscut.st)<-ncpscut.tp
hm <- heatmaply(ncpscut.st, k_row = 7, Colv=NA, hclustfun=xhclust, distfun=xdist, RowSideColors=sidecols, row_dend_left=F, cexRow=0.6, showticklabels = c(F,T), margins=c(1,50,NA,0), label_names=c("ion","time","fracmax"), main=with(sp[pl,],paste("Exp",Experiment," Cuv",Cuvette,":",Species,Population,":",Closed, "closed,",Flrs,"open flowers")))
layout(hm)
}

threshold <- 1e-4
#TODO: time markers

Hybrid

ncps_plot(1)

S. hookeri

ncps_plot(2)
ncps_plot(8)
ncps_plot(9)

S. kaalae

ncps_plot(3)
ncps_plot(4)
ncps_plot(5)
ncps_plot(6)
ncps_plot(7)
trange <- !ncps[[pl]]$pib
prange <- 31:ncol(ncpsp[[pl]])
ncpscut <-  ncpsp[[pl]][trange,prange]
ncpscut <- ncpscut[,colMax(ncpscut)>threshold]
ncpscut[ncpscut<0] <- 0
ncpscut <- log(ncpscut+1)
ncpscut.s <- divmax(ncpscut)
ncpscut.st <- t(as.matrix(ncpscut.s)) 
ncpscut.st3 <- ncpscut.st

ncpscut.tp <- ncps[[pl]]$tp[trange]

ncpscut.frmlas <- frmlas[colnames(ncpscut)]
ncpscut.frmlas <- aggregate(ncpscut.frmlas~round(mass), paste, collapse="/", data=frmla[colnames(ncpscut),])

wmass <- factor(round(frmla[colnames(ncpscut),"mass"]))
ncpscut.st <- as.data.frame(t(ncpscut.s))
ncpscut.st <- aggregate(. ~ wmass, mean, data=ncpscut.st)
rownames(ncpscut.st) <- ncpscut.st$wmass
ncpscut.st$wmass <- NULL

ncpscut.t <- as.data.frame(t(ncpscut))
ncpscut.t <- aggregate(. ~ wmass, mean, data=ncpscut.t)
rownames(ncpscut.t) <- ncpscut.t$wmass
ncpscut.t$wmass <- NULL

colside_max <- range01(rowMax(ncpscut.t)^(1/16))
sidecols <- colorRampPalette(c("white", "purple"))(512)[round(colside_max*511)+1]
xhclust <- function(data) { hclust(data, method="mcquitty") }
xdist <- function(data) { dist(data, method="DTW") }

#heatmap(as.matrix(ncpscut.st), Colv=NA, col=topo.colors(512), hclustfun=xhclust, distfun=xdist, labRow=rownames(ncpscut.st), labCol=NA, RowSideColors=sidecols, margins=c(5,18))
  
#rownames(ncpscut.st)<-ncpscut.frmlas[,2]
colnames(ncpscut.st)<-ncpscut.tp
hm <- heatmaply(ncpscut.st, k_row = 7, Colv=NA, hclustfun=xhclust, distfun=xdist, RowSideColors=sidecols, row_dend_left=F, cexRow=0.6, showticklabels = c(F,T), margins=c(1,50,NA,0), label_names=c("ion","time","fracmax"), main=with(sp[pl,],paste("Exp",Experiment," Cuv",Cuvette,":",Species,Population,":",Closed, "closed,",Flrs,"open flowers")))
layout(hm)
}

image(as.matrix(t(ncpscut.st)))
library(grid)
#TODO: add theshold
divmaxr <- function(x) { sweep(x, 1, apply(x, 1, max)+0.000000001, "/") }
comb <- as.matrix(ncpscut.st1)+as.matrix(ncpscut.st2)+as.matrix(ncpscut.st3)
hord2 <- heatmap(divmaxr(comb[rowSums(comb)>100,]), Colv=NA, cexRow=0.5, hclustfun=xhclust, distfun=xdist)
dev.off()

grid.raster(matrix(rgb(
  red=as.matrix(ncpscut.st3[rowSums(comb)>100,][hord2$rowInd,]),
  green=as.matrix(ncpscut.st2[rowSums(comb)>100,][hord2$rowInd,]), 
  blue=as.matrix(ncpscut.st1[rowSums(comb)>100,][hord2$rowInd,])),
  ncol=dim(ncpscut.st)[2]), interpolate=FALSE)

grid.raster(matrix(rgb(
  red=as.matrix(ncpscut.st3[hord$rowInd,]),
  green=as.matrix(ncpscut.st2[hord$rowInd,]), 
  blue=as.matrix(ncpscut.st1[hord$rowInd,])),
  ncol=dim(ncpscut.st)[2]), interpolate=FALSE)

Discussion

Top volatiles interpretation

Some delay is expected at least the first night if the circadian clock for volatiles was syncronized with ambient sunset. On the second night, the plant may respond to the previous day’s shorter artificial photoperiod.