Introduction to EZFragility package
This is the introduction to how to use the EZFragility package. It
calculates the neural fragility heatmap based on intracranial
electrocorticographic (iEEG) recordings of ictal events. The method is
based on the paper by Li et al. (2017) <doi:
10.23919/ACC.2017.7963378>.
Fragility Method
Key definition of the neural fragility marker
Drug-resistant epilepsy (DRE) is a significant clinical challenge.
Effective treatment requires the synthesis of large amounts of
personalized clinical data with the goal of identifying the region
within a patient’s brain from which seizures originate (the
epileptogenic zone, or EZ). Current standard of care involves
multidisciplinary group meetings and clinician consensus after reviewing
patient data. This lack of a standardized methodology has led to
variable surgical outcomes, measured by seizure freedom/reduction after
surgery. Computational methods to identify the EZ have been developed to
improve outcomes. The neural fragility implemented is the R package is a
promising biomarker of the epileptogenic zone. We provide this R package
so that the method can be tested on extensive data for research.
The hypothesis of the neural fragility is that small changes in
connection strengths at seizure onset zone nodes (SOZ) cause an
imbalance in inhibitory and excitatory connectivity between brain
regions. Then if the SOZ is perturbed, over excitation can occur
manifesting in a seizure. The neural fragility marker propose a metric
to evaluate when connection strength changes cause imbalance and on
which electrodes.To create the fragility marker, a personalized dynamic
model of the brain network from observed iEEG signals is computed in a
discrete number of time windows. Then the neural fragility measures the
degree network nodes are imbalanced by systematically computing the
minimum perturbation required for each electrode’s connections to
produce instability for the entire network.
Load small ictal ieeg data for patient PT01
This example data corresponds to the first seizure of patient PT01
from the Fragility Data Set. pt01EcoG is a matrix with 3000
rows (time points) and 84 columns (electrodes). The acquisition
frequency is 1000Hz and the time around seizure onset is [-1:2]s. The
rows names are the electrodes names. The column names is the time
relative to seizure onset in (s). For more information about the data,
you can use the ?pt01EcoG command.
The patient underwent a surgical treatment and was seizure-free after
the surgery. The electrodes located in the SOZ and the surgically
removed area are marked in sozIndex attribute.
data("pt01EcoG")
## create an epoch object
epoch <- Epoch(pt01EcoG)
epoch
#> 84 x 3000 matrix of doubles:
#>
#> [,1] [,2] [,3] ... [,3000]
#> G1 16650.233 25734.890 37898.546 ... 96182.263
#> G2 -225911.037 -226337.719 -226012.026 ... 108857.517
#> G3 -230090.583 -235515.210 -245789.565 ... 149984.457
#> ... ... ... ... ... ...
#> SLT4 109452.929 103181.323 101214.372 ... 92960.977
#>
#> Time range: -1 to 1.999
#> Use $ to access its methods. see "?`Epoch-method`"
The function visuIEEGData can visualize the iEEG
data.
visuIEEGData(epoch = epoch)
The Epoch object can be subset using the [ operator.
# The electrode names corresponding to the site of the patient's surgery
sozNames <- attr(pt01EcoG, "sozNames")
## Show the electrodes that are marked as SOZ and additional 4 electrodes
display <- c(sozNames, "MLT1", "MLT2", "MLT3", "MLT4")
visuIEEGData(epoch = epoch[display])
## Equivalent to:
## visuIEEGData(epoch = epoch[display, ])
## constrain to the first 100 time points
visuIEEGData(epoch = epoch[display, 1:100])
You can use truncateTime to specify the time range
epochClipped <- truncateTime(epoch, from = -1, to = 0)
visuIEEGData(epoch = epochClipped)
Compute Fragility Matrix
The function calcAdjFrag computes the fragility matrix
for the given iEEG time series. The simplest form is
calcAdjFrag(ieegts, window, step). However, computing
Fragility Matrix can take a bit while. To speed up the computation, we
provided a parallel example below. The parallel feature is built upon
the foreach package. You need to have a parallel backend
registered in your R session beforehand. For example, you can use the
doSNOW package to register a parallel backend. The result
pt01Frag (data("pt01Frag")) from this code snippet
corresponds to the data saved in the package. The function
calcAdjFrag allows the user to tune the parameters window,
step and nSearch (see help ?calcAdjFrag).
## Register a SNOW parallel backend with 4 workers
library(doSNOW)
#> Warning: package 'doSNOW' was built under R version 4.2.3
#> Loading required package: foreach
#> Warning: package 'foreach' was built under R version 4.2.3
#> Loading required package: iterators
#> Warning: package 'iterators' was built under R version 4.2.3
#> Loading required package: snow
#> Warning: package 'snow' was built under R version 4.2.3
cl <- makeCluster(4, type = "SOCK")
registerDoSNOW(cl)
windowNum <- 250
step <- 125
pt01Frag <- calcAdjFrag(epoch = epoch, window = windowNum, step = step, parallel = TRUE, nSearch=100L,progress = TRUE)
# Fragility result
pt01Frag
#>
#> Fragility object
#> -----------------------------------------------------------------
#> Slot | Class | DIM/LEN | Values
#> -----------------------------------------------------------------
#> frag | matrix | [84,23] | 0.6310, 0.6618, 0.7484...
#> R2 | matrix | [84,23] | 0.9703, 0.9848, 0.9936...
#> lambdas | numeric | [23] | 0.0001, 0.0010, 0.0001...
#> startTimes | numeric | [23] | -1.0000, -0.8750, -0.7500...
#> electrodes | character | [84] | G1, G2, G3...
#> -----------------------------------------------------------------
#> Use '$attr' to access the data
# Stop the parallel backend
stopCluster(cl)
Predict the SOZ
High value of the fragility correlate with the SOZ. The function
estimateSOZ estimates the SOZ from the fragility matrix. By
default, it will return the electrode names that are in the top 10% of
the fragility values. In this example,
pt01Frag$startTimes > 0 corresponds to the time windows
that are after the seizure onset. This is not required. However, we
found that this usually improves the results.
soz <- estimateSOZ(pt01Frag[ ,pt01Frag$startTimes > 0])
soz
#> [1] "ATT2" "AD2" "G32" "ATT1" "AD3" "G12" "G13" "SLT3" "PLT5"
Fragility Statistics
The function fragStat computes quantiles, mean and
standard deviation for two electrodes group marked as soz non marked as
soz to check that the high value of fragilities. The significant higher
distribution of the fragility in the electrode group marked as soz can
be used as a prediction of seizure-freedom after resection of these
electrodes. The quantile results can be used in a multi-patient study
with a machine learning model to detect automatically the SOZ Li et
al. (2017) <doi: 10.23919/ACC.2017.7963378>.
stats <- fragStat(pt01Frag, soz)
stats
#>
#> FragStat object (Summary Statistics by Step)
#> ------------------------------------------------------------
#> Slot | Class | DIM/LEN | Values
#> ------------------------------------------------------------
#> qmatrix | matrix | [20,23] | 0.6933, 0.7547, 0.7818...
#> meanSOZ | numeric | [23] | 0.8066, 0.8528, 0.8890...
#> meanRef | numeric | [23] | 0.6738, 0.6780, 0.7161...
#> sdSOZ | numeric | [23] | 0.1203, 0.1024, 0.0622...
#> sdRef | numeric | [23] | 0.1730, 0.1849, 0.1476...
#> ------------------------------------------------------------
#> Use '$attr' to access the data
Plot functions for Fragility Matrix
The function plotFragHeatmap produce a heatmap of the
fragility and allows to visually check the correlation between the SOZ
and sustained high fragility value. Time window around seizure onset
[-1:2]s, which includes the following electrodes:
- Electrodes within the surgery region:
sozNames
- Electrodes outside of the surgery region: “MLT1”, “MLT2”, “MLT3”,
“MLT4”
The parameter sozIndex accepts either electrode names or
their indices in the fragility matrix. If sozIndex is
provided, the function will display the sozIndex electrodes
in blue in the top rows of the heatmap. If sozIndex is not
provided, the function will simply generate a heatmap of the fragility
matrix.
display <- c(sozNames, "MLT1", "MLT2", "MLT3", "MLT4")
plotFragHeatmap(frag = pt01Frag, sozIndex = sozNames)
plotFragDistribution and plotFragQuantile
are similar to the plotFragHeatmap function. They plot the
fragility distribution and quantiles, respectively.
plotFragDistribution(frag = pt01Frag, sozIndex = sozNames)
plotFragQuantile(frag = pt01Frag, sozIndex = sozNames)
