Exponential Regression with Optimal Sample Size Estimates

For System ID on transient signals, stopping criteria for data collection is calculated by referencing precomputed Monte-Carlo simulations.

Exponential regression with Monte-Carlo simulation derived stopping criteria for the optimal sample size

Based On: Optimal design of motor learning experiments informed by Monte-Carlo simulation

Pritesh N. Parmar and James L. Patton

In the field of neuro-rehabilitation, it is often difficult to have patients endure a long session of training, and we seek an estimate of minimum practice trials. In such cases, motor learning measurements across trials are usually exponentially decaying transient signals. Here we employed Monte-Carlo methods to determine the minimum number of samples required from transient responses. We tracked the accuracy of recovery of synthesised data to reveal a prescription for the minimum number of samples for a robust identification of the underlying learning process, given preliminary estimates of the time constant and noise levels. Our results revealed a systematic relationship for the minimum number of samples required from transient signals that can be used as a stopping criteria for data collection. We also evaluated these results by using a past motor learning study and determining the minimum number of required samples (trials) to best estimate learning curves.

_____________________________________________

[p_Nopt,pci_Nopt,stoppingCriteria,best_estimatedParam] = ...

optimalNumSamplesForExpDecayingSignal(tvector,datavector,bootstrapLevel,plotIt)

optimalNumSamplesForExpDecayingSignal calculates stopping criteria for data collection from

transient signals. It uses a sum of exponential functions to estimate transience and then uses

the largest time-constant to estimate the optimal number of samples.

optimalNumSamplesForExpDecayingSignal can handle sparse/missing data. Thus, tvector does not need

to have regular intervals. Missing data can be indicated using NAN within (tvector,datavector)

pair or can simply be omitted.

optimalNumSamplesForExpDecayingSignal uses a bootstrap (resampling or leave-one out/jackknief) method for robust

regression and uses Multistart to find the global optimum.

......................................................

Reference the published article: Parmar, P. N., & Patton, J. L. (2021). Optimal design of motor

learning experiments informed by Monte-Carlo simulation. International Journal of Experimental

Design and Process Optimisation, 6(4), 289-303. <a href="matlab: web('https://doi.org/10.1504/IJEDPO.2021.117248')">Parmar & Patton (2021), IJEDPO</a>.

......................................................

inputs: tvector,datavector = data to be fitted (dimension Nx1 for both; N = number of samples)

bootstrapType = 'resample'(default) or 'jackknief' bootstrap method

plotIt = 1 or 0(default). plots regression fit and stopping criteria estimation.

outputs: p_Nopt = the optimal number of samples (N_opt*) [row 1 = 95% accuracyLevel; row 2 = 99% accuracyLevel]

pci_Nopt = 95% prediction interval for N_opt* [row 1 = 95% accuracyLevel; row 2 = 99% accuracyLevel]

stoppingCriteria = 1 if N_opt* data samples has been collected, and 0 oterwise.

-1 if only non-transient, diverging signals detected.

[row 1 = 95% accuracyLevel; row 2 = 99% accuracyLevel]

best_estimatedParam = fitted parameter estimates; use with generalExponentialFun.m

95% and 99% accuracyLevels denote the regression fit errors at the stoppingCriteria to be within 5% and 1% of their best possible value.

See also: generalExponentialFun.m, expRegressionContinuousMdl3.m, getMLE.m, minTrialFittedMdl_standalone.mat

_________________________________________________________________________

% EXAMPLE#1: first-order

tvector = (1:400)';

datavector = 50*exp(-tvector/60) + 5*randn(size(tvector));

[p_Nopt,pci_Nopt,stoppingCriteria,best_estimatedParam] = optimalNumSamplesForExpDecayingSignal(tvector,datavector,[],1);

% EXAMPLE#2: sparse sampling; second-order

tvector = (1:400)';

datavector = 100*exp(-tvector/60)-100*exp(-tvector/100) + randn(size(tvector));

samplingPercentage = 70/100; % sparse sampling

indToRemove = randsample(1:length(tvector),ceil(length(tvector)*(1-samplingPercentage)),1);

tvector(indToRemove) = nan; datavector(indToRemove) = nan;

[p_Nopt,pci_Nopt,stoppingCriteria,best_estimatedParam] = optimalNumSamplesForExpDecayingSignal(tvector,datavector,[],1);

_________________________________________________________________________

-------------------------------------------------------------------------

Copyright (c) by Dr. Pritesh Parmar, PhD, 2023-07-10.

Changelog:

v1.0 - (07/10/2023) Initial version

v1.1 - (05/09/2025) outputs p_Nopt,pci_Nopt,stoppingCriteria for both 95% and 99% accuracy levels, minor graphical and performance updates

added 'resample' bootstrap method

best bootstrap instance used to address the high sensitivity of residual error

-------------------------------------------------------------------------

capabilities to add in future:

(1) properly detect onset for a exponentially changing signal. Sometimes beginning of a signal may have constant trend before it starts changing exponentially.

(2) expand exponential regression to allow imaginary values.

(3) more general model functions than the sum of exponential models. e.g., y = p1 + (p2+p3*t)*exp(-t/p4) + ... or also sum of sigmoid functions

(4) discrete model function. y_n+1 = p1 + p2*y_n-1 + p3*y_n-2 + ... ; To get time-constant tau for continuous model, use p2 = exp(-(t2-t1)/tau).

(5) N_opt for negative tau; diverging signals