SOCR ≫ TCIU Website ≫ TCIU GitHub ≫

1 Background

The Laplace transform allows us to examine the relations between the space-time and space-kime representations of longitudinal data. The Fourier transformation is a linear operator that maps complex-valued functions of real variables (e.g., space, time domains) to complex valued functions of other real variables (e.g., frequency domain). The Laplace transform is similar. However, it sends complex-valued functions of positive real variables (e.g., time) to complex-valued functions defined on complex variables (e.g., kime).

In this vignette, we introduce the use of inv_kimesurface_transform() and kimesurface_transform() functions in our package. We demonstrate how to do Laplace transform and kimesurface transform with our functions.

require(TCIU)
require(doParallel)
require(cubature)
require(oro.nifti)
require(magrittr)
require(plotly)
require(ggplot2)

2 Laplace Transform (LT) and inverse Laplace Transform (ILT)

2.1 discrete LT and analytical form of LT

Here we first apply the discrete Laplace Transform (LT) function in our package on the sine function with the domain of \([0 : 2\pi]\) to see whether it has the same 2D surface function as the analytic form of LT of sine. The analytic form is \(1/(1+z^2)\).

# For this part of code, we comment it out and import the output plot already generated before to reduce time
# But this part of code can be run successfully. If you are interested, you can try it on your computer!

# discrete Laplace Transform of sine
range_limit = 2
x2 = seq(from = 0, to = range_limit, length.out = 50)[2:50]
  # drop the first row to avoid real part value of 0
y2 = seq(from = 0, to = range_limit, length.out = 50)[2:50]
  # drop the first column to avoid imaginary part value of 0

# Recompute the LT(sin) discretized on lower-res grid
z2_grid = array(dim=c(length(x2), length(y2)))# x2 %o% y2

f_sin <- function(t) { sin(t) }

# kime surface transform
# use parallel computing to speed up code
ncors = 2 # please choose the ncors according to the number of cores your PC has
# it is better that you increase the number of cores used for parallel computing if your computer allows
cl <- makeCluster(ncors)
registerDoParallel(cl)
F = list()
for (i in 1:length(x2) ){
  F[[i]] = 
    foreach(j = 1:length(y2),
            .export='cubintegrate', 
            .packages='cubature') %dopar% {
      TCIU::LT(FUNCT=f_sin, complex(real=x2[i], imaginary = y2[j]))
    }
}
  
stopCluster(cl)
F_vec = lapply(F, unlist)
z2_grid = unlist(do.call(rbind, F_vec))

# explicit form of Laplace Transform of sine
laplace_sine = function(p) { 1/(p^2 + 1) } # Exact Laplace transform of sin(x), continuous function

XY = expand.grid(X=x2, Y=y2)       
complex_xy = mapply(complex, real=XY$X,imaginary=XY$Y)
sine_z =laplace_sine(complex_xy)
dim(sine_z) = c(length(x2), length(y2))# dim(sine_z)  # [1] 49 49
 
# make the two plots in the same plot to compare
lt_ilt_plotly = plot_ly(hoverinfo="none", showscale = FALSE)%>%
  add_trace(z=Re(sine_z)-1, type="surface", surfacecolor=Im(sine_z))  %>%
  add_trace(z = Re(z2_grid), type="surface", opacity=0.7, surfacecolor=Im(z2_grid) )%>%
  layout(title = 
  "Laplace Transform, LT(sin()), Height=Re(LT(sin())), Color=Re(LT(sin())) \n Contrast Exact (Continuous) vs. 
         Approximate (Discrete) Laplace Transform", showlegend = FALSE)

lt_ilt_plotly
sample_save[[1]]

From the plot, we can easily tell that our discrete LT function generate the same 2D surface as the analytic LT function does.

2.2 ILT on discrete LT of sine

Here we first apply the discrete Laplace Transform (LT) function in our package on the sine function with the domain of \([0, 2\pi]\), and then use the inverse Laplace Transform (ILT) function in our package to prove the Laplace Transformation has been converted back to sine.

# For this part of code, we comment it out and import the output plot already generated before to reduce time
# But this part of code can be run successfully. If you are interested, you can try it on your computer!

# discrete Laplace Transform of sine
f_sin = function(t) { sin(t) }
lt_sine = function(z) TCIU::LT(f_sin, z)

# inverse Laplace Transform on the lt_sine
# using parallel computing to speed up code
tvalsn <- seq(0, pi*2, length.out = 20)
cl <- makeCluster(ncors)
registerDoParallel(cl)
sinvalsn <- foreach(t=1:length(tvalsn),
                    .export='cubintegrate', 
                    .packages='cubature')  %dopar% { 
  TCIU::ILT(FUNCT=lt_sine, t=tvalsn[t])
  }
stopCluster(cl)
sinvalsn = unlist(sinvalsn)

# make the plot of the result from ILT 
# to see whether it still looks like sine
sinvalsn_df2 <- as.data.frame(cbind(Re=Re(sinvalsn),Im=Im(sinvalsn),
                                    Sin=sin(tvalsn), time_points=tvalsn))
lt_ilt_sine = ggplot(sinvalsn_df2, aes(x=time_points))+
  geom_line(aes(y=Re, color="Real"), linetype=1, lwd=2) +
  geom_line(aes(y = Sin, color="Sin"), linetype=2, lwd=1) + 
  scale_color_manual(name="Index",
                     values = c("Real"="steelblue", "Sin"="darkred"))+
  labs(title = "Original fMRI Time-series f(t)=sin(t) and \n Reconstructed f'(t)=ILT(F)=ILT(discrete LT(f))",
       subtitle = bquote("F" ~ "=" ~ "discrete LT(sine)")) + 
  xlab("Time") + ylab("fMRI Image Intensities (f and f')") +
  theme_grey(base_size = 16) +
  theme(legend.title = element_text(size=14, color = "black", face="bold"),
        panel.grid.minor.y = element_blank(),
        panel.grid.major.y = element_blank(),
        plot.title = element_text(hjust = 0.5),
        plot.subtitle = element_text(hjust = 0.5))

lt_ilt_sine
sample_save[[2]]