Location: BG_Ca_leak @ 9e42555de571 / matlab_parameter_fitting / K_kappa_ABmassAction.m

Author:

Shelley Fong <s.fong@auckland.ac.nz>

Date:

2022-03-03 13:59:42+13:00

Desc:

Changing Xi to X_i, and adding python param finder code

Permanent Source URI:

https://models.physiomeproject.org/workspace/835/rawfile/9e42555de571e0ea1b4cad67686a66dfba8cfcd8/matlab_parameter_fitting/K_kappa_ABmassAction.m

% finding K, kappa for simple mass-action reaction of Ca leaking from SR
% into cytosol.
% Not gated, and not voltage dependent.


W_i = 25.8;
W_SR = W_i * 0.088235;
z = 1; % for Ca, z is 2  
N_A = 6.022e23;
x_Kp_channel = 725/N_A*1e15; % unit fmol
x_Nab_channel = x_Kp_channel; % number of channels in myocyte [=] fmol
F = 96485;


%% stoichiometric matrices
% include kf and kr in matrices
N_f = [1; 0]; % [Ca_SR, Cai]
N_r = [0; 1];

N_fT = transpose(N_f);
N_rT = transpose(N_r);
N = N_r - N_f;

num_cols = size(N,2); % number of reactions
I = eye(num_cols); 

M = [I N_fT; I N_rT];

% k+ = k- = K_leak (on ClancyRudy original model)
K_leak = 0.5833333;
kf = [K_leak]; 
kr = [K_leak];

k = [kf;kr];
W = [1;W_SR;W_i];

% solve lambda
lambdaW = exp(pinv(M)*log(k));
lambda = lambdaW./W;
kappa = lambda(1:size(N,2));
K = lambda(size(N,2)+1:end);

save(['K_kappa.mat'],'kappa','K');

% Checks for kappa and K
N_rref = rref(N);
R_mat = null(N,'r');

K_eq = kf./kr;
zero_est = R_mat'*K_eq;

k_est = exp(M*log(lambdaW));
diff = sum(abs((k-k_est)./k));

%% display text for cellml
rname = {'leak'};
Kname = {'CaSR','Cai'};
for ik = 1:length(kappa)
    fprintf('var kappa_%s: fmol_per_sec {init: %g, pub: out};\n', rname{ik}, kappa(ik));
end
for ik = 1:length(K)
    fprintf('var K_%s: per_fmol {init: %g, pub: out};\n', Kname{ik}, K(ik));
end
if true
    % component import parameters
    for ik = 1:length(kappa)
        fprintf('var kappa_%s: fmol_per_sec {pub: in};\n', rname{ik});
    end
    for ik = 1:length(K)
        fprintf('var K_%s: per_fmol {pub: in};\n', Kname{ik});
    end
    
    % mapping
    for ik = 1:length(kappa)
        fprintf('vars kappa_%s and kappa_%s;\n', rname{ik}, rname{ik});
    end
    for ik = 1:length(K)
        fprintf('vars K_%s and K_%s;\n', Kname{ik}, Kname{ik});
    end    
end