SOCR ≫ TCIU Website ≫ TCIU GitHub ≫

This Spacekime TCIU Learning Module presents the core elements of spacekime analytics including:

  • Import of repeated measurement longitudinal data,
  • Numeric (stitching) and analytic (Laplace) kimesurface reconstruction from time-series data,
  • Forward prediction modeling extrapolating the process behavior beyond the observed time-span \([0,T]\),
  • Group comparison discrimination between cohorts based on the structure and properties of their corresponding kimesurfaces. For instance, statistically quantify the differences between two or more groups,
  • Unsupervised clustering and classification of individuals, traits, and other latent characteristics of cases included in the study,
  • Construct low-dimensional visual representations of large repeated measurement data across multiple individuals as pooled kimesurfaces (parameterized 2D manifolds),
  • Statistical comparison, topological quantification, and analytical inference using kimesurface representations of repeated-measurement longitudinal data.

1 Preliminary setup

TCIU and other R package dependencies …

library(TCIU)
library(DT)
library("R.matlab")
# library(AnalyzeFMRI) # https://cran.r-project.org/web/packages/AnalyzeFMRI/index.html
library(ggfortify)
library(plotly)
library(manipulateWidget)
library(transport)
library(shapes)

2 Longitudinal Data Import

In this case, we are just loading some exemplary fMRI data, which is available here.

# pathname <- file.path("./test_data", "subj1_run1_complex_all.mat")
mat1 = readMat("./test_data/subj1_run1_complex_all.mat")
bigim1 = mat1$bigim[,64:1,,]
bigim1_mod = Mod(bigim1) # Modulus
smoothmod = GaussSmoothArray(bigim1_mod,sigma = diag(3,3))
# dim(bigim1) = 64 64 40
# bigim1 contains the complex image space
# dimensions are 64x*64y*40z*160t, corresponding to x,y,z,time
load("./test_data/mask.rda") # load the 3D nifti data of the mask
load("./test_data/mask_label.rda") # load the 3D nifti data of the mask with labels
load("./test_data/mask_dict.rda") # load the label index and label name
label_index = mask_dict$index
label_name = as.character(mask_dict$name)
label_mask = mask_label
load("./test_data/hemody.rda") # load the hemodynamic contour of the brain

Alternatively, we can simulate synthetic ON/OFF fMRI data that can be saved as a pair of arrays ON (stimulus) and OFF (rest) conditions.

library(plotly)
library(dplyr)

# 1. Simulate On/Off fMRI time-series data
simulate_fmri_data <- function(n_time = 80, n_runs = 5, subjects = 20) {
  set.seed(123)
  
  # Example: 'ON' data
  arr_on <- array(
    rnorm(n_time * n_runs * subjects, mean=0, sd=1),
    dim = c(n_time, subjects, n_runs)
  )
  # 'OFF' could be something different. For demonstration, let's shift it
  arr_off <- array(
    rnorm(n_time * n_runs * subjects, mean=1, sd=1),
    dim = c(n_time, subjects, n_runs)
  )
  
  list(ON = arr_on, OFF = arr_off)
}

# 2. Display the On and Off series
# This function converts the 3D array data into a tidy data frame format
# Calculates the mean time series for each condition (ON/OFF)
# Creates an interactive plot with:
# 
# Individual run data shown with low opacity
# Mean condition data shown with thicker lines
# Interactive tooltips with details on hover
# Legend distinguishing between conditions and runs
# Appropriate axis labels and title
# Function to create interactive plotly visualization of fMRI data
plot_fmri_timeseries <- function(fmri_data, subjects = c(1), title="On/Off fMRI Time-Series Data") {
  # Extract data from the simulation
  on_data <- fmri_data$ON[ , subjects, ]
  off_data <- fmri_data$OFF[ , subjects, ]
  
  # fix the dimensions
  on_data  <- array(on_data, dim = c(dim(on_data)[1], length(subjects), dim(on_data)[2]))
  off_data <- array(off_data, dim = c(dim(off_data)[1], length(subjects), dim(off_data)[2]))
    
  n_time <- dim(on_data)[1]
  n_runs <- dim(on_data)[3]
  
  # Create a data frame for plotting
  plot_data <- data.frame()
  
  # Process ON data
  for (run in 1:n_runs) {
    temp_df <- data.frame(Time = 1:n_time, Value = on_data[, 1, run],
                          Condition = "ON", Run = paste("Run", run))
    plot_data <- rbind(plot_data, temp_df)
  }
  
  # Process OFF data
  for (run in 1:n_runs) {
    temp_df <- data.frame(Time = 1:n_time, Value = off_data[, 1, run],
                          Condition = "OFF", Run = paste("Run", run))
    plot_data <- rbind(plot_data, temp_df)
  }
  
  # Calculate mean values for each condition across runs
  mean_data <- plot_data %>%
    group_by(Time, Condition) %>%
    summarize(Mean_Value = mean(Value), .groups = 'drop')
  
  # Create interactive plot
  p <- plot_ly() %>%
    # Add individual run traces with lower opacity
    add_trace(data = plot_data, x = ~Time, y = ~Value, color = ~Condition,
               type = 'scatter', mode = 'lines', linetype = ~Run, opacity = 0.3,
               showlegend = TRUE, hoverinfo = 'text', text = ~paste('Condition:',
                  Condition, '<br>Run:', Run, '<br>Time:', Time, '<br>Value:', round(Value, 3))
    ) %>%
    # Add mean traces with higher opacity and wider lines
    add_trace(data = mean_data, x = ~Time, y = ~Mean_Value, color = ~Condition,
               type = 'scatter', mode = 'lines', line = list(width = 4),
               name = ~paste(Condition, "Mean"), hoverinfo = 'text',
               text = ~paste('Condition:', Condition, 'Mean', 
                   '<br>Time:', Time, '<br>Value:', round(Mean_Value, 3))
    ) %>%
    layout(title = title,
      xaxis = list(title = "Time Point"), yaxis = list(title = "Signal Value"),
      legend = list(orientation = 'h', y=-0.1), hovermode = "closest")
  return(p)
}

# 2.1 Generate simulated data
fmri_data <- simulate_fmri_data(n_time = 80, n_runs = 5, subjects = 20)
# 2.2 Create interactive plot
plot_fmri_timeseries(fmri_data=fmri_data, subjects = c(1), title="Simulated On/Off fMRI Time-Series Data")


# 3. The triple Kimesurface visualization function
visualize_kimesurfaces_tripleREV <- function(fmri_data) {
  t_min <- 1;     t_max <- dim(fmri_data$ON)[1]
  th_min<- 1;     th_max<- dim(fmri_data$ON)[2]*dim(fmri_data$ON)[3]
  
  t_seq     <- c(t_min:t_max)
  theta_seq <- c(th_min:th_max)
  
  # flatten the repeats within subject over runs (5) and across subjects (20) == th_max
  # into a single dimensions, random sampling for theta,
  reshapedOn <- matrix(fmri_data$ON, nrow = th_max, ncol = t_max)
  reshapedOff <- matrix(fmri_data$OFF, nrow = th_max, ncol = t_max)
  
  # 2D surface in (t,theta) space
  fig <- plot_ly() %>%
    # First: ON
    add_surface(
      z = reshapedOn,
      colorscale = "Viridis",
      showscale  = FALSE,
      name       = "ON",
      scene      = "scene1"
    ) %>%
    layout(
      scene = list(
        domain     = list(x = c(0, 0.33), y = c(0, 1)),
        aspectmode = "cube",
        xaxis = list(title = "time (t)"),
        yaxis = list(title = "Runs * Participants (θ)"),
        zaxis = list(title = "Intensity")
      )
    ) %>%
    # Second: OFF
    add_surface(
      z = reshapedOff,
      colorscale = "Viridis",
      showscale  = FALSE,
      name       = "OFF",
      scene      = "scene2"
    ) %>%
    layout(
      scene2 = list(
        domain     = list(x = c(0.33, 0.66), y = c(0, 1)),
        aspectmode = "cube",
        xaxis = list(title = "time (t)"),
        yaxis = list(title = "Runs * Participants (θ)"),
        zaxis = list(title = "Intensity")
      )
    ) %>%
    # Third: Overlay
    add_surface(
      z = reshapedOn,
      colorscale = "Viridis",
      opacity    = 0.5,
      showscale  = FALSE,
      name       = "ON (overlay)",
      scene      = "scene3"
    ) %>%
    add_surface(
      z = reshapedOff,
      colorscale = "Portland",
      opacity    = 0.8,
      showscale  = FALSE,
      name       = "OFF (overlay)",
      scene      = "scene3"
    ) %>%
    layout(
      scene3 = list(
        domain     = list(x = c(0.66, 1), y = c(0, 1)),
        aspectmode = "cube",
        xaxis = list(title = "time (t)"),
        yaxis = list(title = "Runs * Participants (θ)"),
        zaxis = list(title = "Intensity")
      )
    ) %>%
    layout(
      title      = "Kimesurface Visualization of fMRI ON and OFF",
      showlegend = TRUE,
      legend     = list(x=0.85, y=0.1),
      annotations= list(
        list(showarrow=FALSE, text='(ON)',  xref='paper', yref='paper',
             x=0.15, y=0.95, xanchor='center', yanchor='bottom',
             font=list(size=16)),
        list(showarrow=FALSE, text='(OFF)', xref='paper', yref='paper',
             x=0.5,  y=0.95, xanchor='center', yanchor='bottom',
             font=list(size=16)),
        list(showarrow=FALSE, text='(OVERLAY)',
             xref='paper', yref='paper', x=0.85, y=0.95,
             xanchor='center', yanchor='bottom',
             font=list(size=16))
      )
    ) %>%
    # Sync the cameras
    onRender("
      function(el) {
        el.on('plotly_relayout', function(d) {
          const camera = Object.keys(d).filter((key) => /\\.camera$/.test(key));
          if (camera.length) {
            const scenes = Object.keys(el.layout).filter((key) => /^scene\\d*/.test(key));
            const new_layout = {};
            scenes.forEach(key => {
              new_layout[key] = {...el.layout[key], camera: {...d[camera]}};
            });
            Plotly.relayout(el, new_layout);
          }
        });
      }
    ")
  
  fig
}

# 4. Display the kime-surfaces
fmri_data <- simulate_fmri_data()
fig <- visualize_kimesurfaces_tripleREV(fmri_data)
fig

3 Time-series graphs

3.1 Interactive time-series visualization

The TCIU function fmri_time_series() is used to create four interactive time series graphs for the real, imaginary, magnitude, and phase parts for the fMRI spacetime data. We can also add a reference plotly object to the plot. This function is based on the GTSplot function from package TSplotly.

3.1.1 Example fMRI(x=4, y=42, z=33, t)

# reference_plot = sample[[4]]
fmri_time_series(bigim1, c(44,42,33), is.4d = TRUE) #, ref = reference_plot)

4 Kime-series/kimesurfaces (spacekime analytics protocol)

The following examples connect to several ongoing spacekime analytics, including kimesurface analyticity, kimesurface regularization, and distribution-time dynamics.

The first example uses a time-series simulation to illustrate how to transform the fMRI time-series at a fixed voxel location into a kimesurface (kime-series).

Notes:

4.1 Pseudo-code

  • Randomly generate \(8\) \(phi=\phi\) kime-phases for each of the \(10\) time radii. This yields an \(8\times 10\) array (phi_8_vec) of kime phase directions. These directions can be obtained by different strategies, e.g., (1) uniform or Laplace distribution sampling over the interval \([-\pi:\pi)\), (2) randomly sampling with/without replacement from known kime-phases obtained from similar processes, etc.
  • Optionally, order all kime-phases (rows) from small to large for each column.
  • Fix the \(\nu=(x,y,z)\) voxel location and extract the fMRI time-course \(fMRI_{\nu}(t), \forall 1\leq t\leq 160\).
  • For binary stimuli (e.g., activation (ON) and rest (OFF) event-related design), we can split the 160-long fMRI series into \(80\) ON (Activation) and \(80\) OFF (rest) states, or sub-series.
  • Construct a data-frame with \(160\) rows and \(4\) columns; time (\(1:10\)), phases (\(8\)), states (\(2\)), and fMRI_value (Complex or Real intensity).
  • Convert the long DF representing fMRI_ON and fMRI_OFF from their native (old) polar coordinates to the (new) Cartesian coordinates, using polar transformations.
  • Check for visual (graphical) and numeric differences between the fMRI intensities during the ON vs. OFF states
  • Spatially smooth (blur) the matrices (in 2D) to reduce noise make them more representative of the process. May also need to temper the magnitude of the raw fMRI intensities, which can have a large range.
  • Generate the plot_ly text labels that will be shown on mouse hover (pop-up dialogues) over each kimesurface/kime-series. These text-labels are stored in Cartesian coordinates \((-10\leq x\leq 10,-10\leq y\leq 10)\), but are computed using the polar coordinates \((1\leq t\leq 10,-\pi\leq \phi<\pi)\) and the polar-to-Cartesian transformation. The labels are \(21\times21\) arrays including the triple \((t, \phi, fMRIvalue)\). Separate text-labels are generated for each kimesurface (ON vs. OFF stimuli).
  • Generate the \(21\times21\) kime-domain Cartesian grid by polar-transforming the polar coordinates \((t,\phi)\) into Cartesian counterparts \((x,y)\).
  • Interpolate the fMRI intensities (natively anchored at \((t,\phi)\)) to \(fMRI(x,y), \forall -11\leq x,y \leq 11, \sqrt{x^2+y^2}\leq 10\).
  • Use plot_ly to display in 3D the kime-series as 2D manifolds (kimesurfaces) over the Cartesian domain.

4.2 Function main step: Time-series to kimesurfaces Mapping

4.2.1 Generate the kime-phases

Randomly generate \(8\) \(phi=\phi\) kime-phases for each of the \(10\) time radii. This yields an \(8\times 10\) array (phi_8_vec) of kime phase directions. These directions can be obtained by different strategies, e.g., (1) uniform or Laplace distribution sampling over the interval \([-\pi:\pi)\); (2) randomly sampling with/without replacement from known kime-phases obtained from similar processes; (3) use an \(AR(p)\) autoregressive model with memory (affected by prior samples at earlier \(p\) times), etc.

Optionally, order all kime-phases (rows) from small to large for each column.

# plot Laplacian
# ggfortify::ggdistribution(extraDistr::dlaplace, seq(-pi, pi, 0.01), m=0, s=0.5)

t <- seq(-pi, pi, 0.01)
f_t <- extraDistr::dlaplace(t, mu=0, sigma=0.5)
plot_ly(x=t, y=f_t, type="scatter", mode="lines", name="Laplace Distribution") %>%
  layout(title="Laplace Distribution: Kime-Phase Prior")
# randomly generate 8 phi kime-phases for each of the 10 time radii
phi_8_vec <- matrix(NA, ncol=10, nrow = 8)
for (t in 1:10) { 
  # for a given t, generate 8 new phases
  set.seed(t);
  phi_8_vec[ ,t] <-
    extraDistr::rlaplace(8, mu=0, sigma=0.5)
  # rank-order the phases for consistency
  # within the same foliation leaf
  phi_8_vec[ ,t] <- sort(phi_8_vec[ ,t])
  # force phases in [-pi: pi)
  for (i in 1:8) {
    if (phi_8_vec[i,t] < -pi) 
      phi_8_vec[i,t] <- -pi
    if (phi_8_vec[i,t] >= pi) 
      phi_8_vec[i,t] <- pi
  }
}
# order all kime-phases (rows) from small to large for each column
# phi_8_vec_ordered <- apply(phi_8_vec, 2, sort)

4.2.2 Structural Data Preprocessing

Here should be included any study-specific data preprocessing. In the case of fMRI series, we may need to split off the individual repeated measurement fMRI time-series from the master (single) fMRI time-series according to the specific event-related fMRI design.

For simplicity, consider a simulated binary stimulus paradigm (e.g., activation (ON) and rest (OFF) event-related design). We can split the \(160\)-timepoint fMRI series into \(80\) ON (Activation) and \(80\) OFF (rest) states, or sub-series, consisting of \(8\) repeats, each of length \(10\) time points, where each time point corresponds to about \(2.5\ s\) of clock time.

fMRI_ON<-bigim1_mod[40,42,33,][c(rep(TRUE,10),rep(FALSE,10))]
fMRI_OFF<-bigim1_mod[40,42,33,][c(rep(FALSE,10),rep(TRUE,10))]

4.2.3 Intensity Data Preprocessing

In practice, some spatial smoothing (blurring) the 1D time-series or their 2D array (tensor representations as 2D matrices (\(row=repeat\), \(column=time\)) to reduce noise and make the data more natural (low-pass filtering, avoiding high-pitch noise). Sometimes, we may also need to temper the magnitude of the raw time-series intensities, which can have a large range.

# Convert the long DF representing fMRI_ON and fMRI_OFF from polar coordinates to Cartesian coordinates
library(spatstat)

matrix_ON <- matrix(0, nrow = 21, ncol = 21) 
matrix_OFF <- matrix(0, nrow = 21, ncol = 21) 
for (t in 1:10) {
  for (p in 1:8) {
    x = 11+t*cos(phi_8_vec[p,t])
    y = 11+t*sin(phi_8_vec[p,t])
    matrix_ON[x,y]  <- fMRI_ON[(p-1)*10 +t] # What is matrix_ON?
    matrix_OFF[x,y] <- fMRI_OFF[(p-1)*10 +t]
  }
}
# smooth/blur the matrices
# matrix_ON_smooth <- (1/10000)*as.matrix(blur(as.im(matrix_ON), sigma=0.5))
# matrix_OFF_smooth <- (1/10000)*as.matrix(blur(as.im(matrix_OFF), sigma=0.5))
matrix_ON_smooth <- (1/10000)*as.matrix(blur(as.im(matrix_ON), sigma=0.5))
matrix_OFF_smooth <- (1/10000)*as.matrix(blur(as.im(matrix_OFF), sigma=0.5))

4.2.4 Generate plotly labels

Generate the plot_ly text labels that will be shown over the graph, upon mouse hovering (pop-up dialogues) over each kimesurface/kime-series. These text-labels are stored in Cartesian coordinates \((-10\leq x\leq 10,-10\leq y\leq 10)\), but are computed using the polar coordinates \((1\leq t\leq 10,-\pi\leq \phi<\pi)\) and the polar-to-Cartesian transformation. The labels are \(21\times 21\) arrays including the triple \((t, \phi, fMRIvalue)\). Separate text-labels are generated for each kimesurface (ON vs. OFF stimuli).

# fix the plot_ly Text Labels
x <- vector()
y <- vector()
i <- 1
for (t in 1:10) {
  for (p in 1:8) {
    x[i] = 11+t*cos(phi_8_vec[p,t])
    y[i] = 11+t*sin(phi_8_vec[p,t])
    i <- i+1
  }
}
# hoverText <- cbind(x=1:21, y=1:21, height=as.vector(t(matrix_ON_smooth))) # tail(mytext)
# custom_txt <- matrix(NA, nrow=21, ncol=21)
# hoverTextOFF <- cbind(x=1:21, y=1:21, height=as.vector(t(matrix_OFF_smooth))) # tail(mytext)
# custom_txtOFF <- matrix(NA, nrow=21, ncol=21)
# 
# for (x in 1:21) {
#    for (y in 1:21) {
#      t = sqrt((x-11)^2 + (y-11)^2)
#      p = atan2(y-11, x-11)
#      custom_txt[x,y] <- paste(' fMRI: ', round(hoverText[(x-1)*21+y, 3], 3),
#                     '\n time: ', round(t, 0),
#                     '\n phi: ', round(p, 2)) 
#      custom_txtOFF[x,y] <- paste(' fMRI: ', round(hoverTextOFF[(x-1)*21+y, 3], 3),
#                     '\n time: ', round(t, 0),
#                     '\n phi: ', round(p, 2)) 
#    }
# }

4.2.5 Cartesian space interpolation

Interpolate the fMRI intensities, natively anchored at polar (kime) coordinates) \(\left (\underbrace{t}_{time},\underbrace{\phi}_{repeat} \right)\), into Cartesian coordinates \(fMRI(x,y), \forall -11\leq x,y \leq 11, x^2+y^2\leq 10\). Note that for specificity, we hard-coded polar-grid parameterization to time \(t\in\{1,2,3, \cdots,10\},\ |T|=10\) and phase \(\phi=seq\left (from=-\pi, to=\pi, step=\frac{2\pi}{20}\right )\in\{-3.1415927, -2.8274334, \cdots, 0.0,\cdots, 2.8274334, 3.1415927 \},\ |\Phi|=21\).

# xx2 <- 11 + c(-10:10) %o% cos(seq(-pi, pi, 2*pi/20))
# yy2 <- 11 + c(-10:10)  %o% sin(seq(-pi, pi, 2*pi/20))
xx2 <- 11 + seq(0,10,1/2) %o% cos(seq(-pi, pi, 2*pi/20))
yy2 <- 11 + seq(0,10,1/2)  %o% sin(seq(-pi, pi, 2*pi/20))
#zz2 <- as.vector(matrix_ON_smooth) %o% rep(1, 21*21)
zz2 <- matrix_ON_smooth
ww2 <- matrix_OFF_smooth
dd2 <- matrix_ON_smooth-matrix_OFF_smooth

#plot 2D into 3D and make the text of the diameter (time), height (r), and phase (phi)
f <- list(family = "Courier New, monospace", size = 18, color = "black")
x <- list(title = "k1", titlefont = f)
y <- list(title = "k2", titlefont = f)
z <- list(title = "fMRI Kime-series", titlefont = f)
# matrix_ON_smooth,matrix_OFF_smooth is indexed in Cartesian coordinates
# Reinterpolate the fMRI intensities according to Polar grid
library(akima)
ON_transformed <- matrix(0, nrow = 21, ncol = 21)
OFF_transformed <- matrix(0, nrow = 21, ncol = 21)
ON_OFF_transformed <- matrix(0, nrow = 21, ncol = 21)

cart_x <- as.vector(rep(seq(1,21),21))
cart_y <- as.vector(sort(rep(seq(1,21),21)))
cart_z_on <- as.vector(matrix_ON_smooth)
cart_z_diff <- as.vector(matrix_ON_smooth-matrix_OFF_smooth)
cart_z_off <- as.vector(matrix_OFF_smooth)

for(i in 1:21){
  for(j in 1:21){
    #interpolations
    int_res_on <- interp(cart_x,cart_y,cart_z_on,xx2[i,j],yy2[i,j])
    int_res_off <- interp(cart_x,cart_y,cart_z_off,xx2[i,j],yy2[i,j])
    int_res_diff <- interp(cart_x,cart_y,cart_z_diff,xx2[i,j],yy2[i,j])
    #insert data
    ON_transformed[i,j] = as.numeric(int_res_on$z)
    OFF_transformed[i,j] = as.numeric(int_res_off$z)
    ON_OFF_transformed[i,j] = as.numeric(int_res_diff$z)
    # if None set to 0
    if(is.na(ON_transformed[i,j])){
      ON_transformed[i,j] = 0
    }
    if(is.na(ON_OFF_transformed[i,j])){
      ON_OFF_transformed[i,j] = 0
    }
    if(is.na(OFF_transformed[i,j])){
      ON_OFF_transformed[i,j] = 0
    }
  }
}
# plot_ly(x = xx2,y = yy2,z = ON_transformed) %>% add_surface()

4.2.6 Cartesian representation

Generate the \(21\times 21\) kime-domain Cartesian grid hovering texts.

# hoverText <- cbind(x=1:21, y=1:21, height=as.vector(t(matrix_ON_smooth))) # tail(mytext)
custom_txt <- matrix(NA, nrow=21, ncol=21)
# hoverTextOFF <- cbind(x=1:21, y=1:21, height=as.vector(t(matrix_OFF_smooth))) # tail(mytext)
custom_txtOFF <- matrix(NA, nrow=21, ncol=21)
custom_txt_DIFF <- matrix(NA, nrow=21, ncol=21)
custom_txt_ON_OFF <- matrix(NA, nrow=21, ncol=21)

for (xdir in 1:21) {
   for (ydir in 1:21) {
     custom_txt[xdir,ydir] <- paste(' fMRI: ', round(ON_transformed[xdir,ydir], 3),
                    '\n time: ', round((xdir-1)/2, 0),
                    '\n phi: ', round(-pi+pi/10*(ydir-1), 2))
     custom_txtOFF[xdir,ydir] <- paste(' fMRI: ', round(OFF_transformed[xdir,ydir], 3),
                    '\n time: ', round((xdir-1)/2, 0),
                    '\n phi: ', round(-pi+pi/10*(ydir-1), 2))
     custom_txt_DIFF[xdir,ydir] <- paste(' fMRI: ', round(ON_transformed[xdir,ydir]-OFF_transformed[xdir,ydir], 3),
                    '\n time: ', round((xdir-1)/2, 0),
                    '\n phi: ', round(-pi+pi/10*(ydir-1), 2))
     custom_txt_ON_OFF[xdir,ydir] <- paste(' fMRI: ', round(ON_OFF_transformed[xdir,ydir], 3),
                    '\n time: ', round((xdir-1)/2, 0),
                    '\n phi: ', round(-pi+pi/10*(ydir-1), 2))
   }
}
plot_ly(x = xx2,y = yy2,z = ON_transformed, type = "surface", colors=c("#FFFFFF","#0000FF"),
          text = custom_txt, hoverinfo = "text", showlegend = FALSE) %>%
    add_trace(x=11, y=11, z=0:0.15, type="scatter3d", mode="lines",
              line = list(width = 10, color="red"), name="Space(x)",
              hoverinfo="none", showlegend = FALSE) %>%
    layout(dragmode = "turntable", title = "ON kimesurface/Kime-Series at a fixed voxel location",scene = list(xaxis = x, yaxis = y, zaxis = z), showlegend = FALSE)
    # 
plot_ly(x = ~xx2, y = ~yy2, z = ~OFF_transformed, type = "surface",colors=c("#FFFFFF","#0000FF"),   # scatterpolar
          text = custom_txtOFF, hoverinfo = "text", showlegend = FALSE) %>% 
    #add_trace(x=~xx2, y=~yy2, z=~ww2, colors = c("blue", "yellow"),
    #          type="surface", text = custom_txtOFF, hoverinfo = "text",
    #          opacity=0.3, showscale = FALSE, showlegend = FALSE) %>%
    # trace the main Z-axis
    add_trace(x=11, y=11, z=0:0.15, type="scatter3d", mode="lines", 
              line = list(width = 10, color="red"), name="Space(x)", 
              hoverinfo="none", showlegend = FALSE) %>%
    layout(dragmode = "turntable", title = "OFF kimesurface/Kime-Series at a fixed voxel location",
           scene = list(xaxis = x, yaxis = y, zaxis = z), showlegend = FALSE)