You are an expert in writing code in many different languages. Your goal is to perform code completion for the following code, keeping in mind the rest of the code and the file meta data.
    Metadata details:
    filename: DEM_I3_and_TS.m
    path of file: ./repos/bspm/thirdparty/spm12/toolbox/DEM
    
    the code of the file until where you have to start completion: function DEM_I3_and_TS

    
function DEM_I3_and_TS
%--------------------------------------------------------------------------
% This routine is a phenomenological examination of the relationship
% between conditional mutual information (i.e.,expected intrinsic mutual
% information and the exponential divergence of trajectories as the
% Rayleigh parameter of a Lorenz attractoris increased (through a pitchfork
% bifurcation and subsequent (subcritical) Hopf bifurcation. The
% (stochastic) Lorentz system is integrated for different values of the
% relay (control) parameter. The nonequilibrium steady-state density is
% then estimated by bidding into a discrete state space; while the
% bifurcations are characterised in terms of the maximal Lyapunov exponent.
% The key thing to observe is the increase in  conditional mutual
% information following the Hopf bifurcation and implicit exponential
% divergence of trajectories. This is scored by the maximal Lyapunov
% exponent crossing zero.


% generative model
%==========================================================================                       % switch for demo
spm_figure('GetWin','DEM'); clf

% flow and Jacobian functions
%--------------------------------------------------------------------------
f  = @(x,v,P,G) v(:) + [-P(1) P(1) 0; P(3) -1 -x(1); x(2) 0 P(2)]*x/64;
Df = @(x,v,P,G) [-P(1) P(1) 0; P(3) -1 -x(1); x(2) 0 P(2)]/64 + ...
                [0 0 0; -x(3) 0 0; 0 x(1) 0]/64;

% model with initial states
%--------------------------------------------------------------------------
G.x  = [1; 1; 24];
G.f  = f;

% set up
%--------------------------------------------------------------------------
b     = 32;                         % number of bins for density estimation
T     = 2^18;                       % length of trajectory
W     = 32;                        % precision of intrinsic fluctuations
P     = exp((-16:32)*log(32)/32);   % Rayleigh parameter range
for k = 1:length(P)
    
    % integrated timeseries
    %----------------------------------------------------------------------
    U.u   = randn(T + 256,3)/W;
    Pk    = [10; -8/3; P(k)];
    t     = spm_int_L(Pk,G,U);
    
    % remove intial transients
    %----------------------------------------------------------------------
    t     = t(256:end