Skip to content

Commit ca96f67

Browse files
alexpopingaalexpopinga
committed
Adapt Dusp1 model for GR-beta
1 parent b511e89 commit ca96f67

File tree

1 file changed

+64
-51
lines changed

1 file changed

+64
-51
lines changed

‎WorkSpace/EricModel/GR_2_Electric_Boogaloo.m

+64-51
Original file line numberDiff line numberDiff line change
@@ -368,8 +368,10 @@
368368
ModelGR_b.stoichiometry = [-1,1,-1,1,0;...
369369
1,0,0,-1,-1];
370370
ModelGR_a.parameters = ({'kn2c',0.01;'dc',1e-5;'dn',1e-6;'ba1',14e-5;...
371-
'MDex',5;'gDex',0.003;'kcn0',0.005;'kcn1',0.02;'Dex0',100});
372-
ModelGR_b.parameters = ({'kn2c',0.01;'dc',1e-5;'dn',1e-6;'kb1',0.01;'bb1',14e-5});
371+
'MDex',5;'gDex',0.003;'kcn0',0.005;'kcn1',0.02;...
372+
'koff',0.1;'kon',0.1;'kr',1;'dr',0.02;'Dex0',100});
373+
ModelGR_b.parameters = ({'kn2c',0.01;'dc',1e-5;'dn',1e-6;'kb1',0.01;'bb1',14e-5;...
374+
'koff',0.1;'kon',0.1;'kr',1;'dr',0.02});
373375
ModelGR_a.inputExpressions = {'IDex','Dex0*exp(-gDex*t)'};
374376
disp('The GR-alpha model is: ')
375377
ModelGR_a.summarizeModel
@@ -491,8 +493,8 @@
491493
ModelGRfit_a{i} = ModelGR_a.loadData("../EricModel/EricData/Gated_dataframe_Ron_020224_NormalizedGR_bins.csv",...
492494
{'nucGR_a','normgrnuc';'cytGR_a','normgrcyt'},...
493495
{'Dex_Conc',GRfitCases{i,2}});
494-
ModelGRfit_a{i}.parameters(9,:) = {'Dex0', str2num(GRfitCases{i,1})};
495-
ModelGRparameterMap_a(i) = {(1:8)}; % Parameters 1:8 refers to all ModelGR_a parameters with the exception of Dex0.
496+
ModelGRfit_a{i}.parameters(13,:) = {'Dex0', str2num(GRfitCases{i,1})};
497+
ModelGRparameterMap_a(i) = {(1:8)};
496498
end
497499

498500
% Make Guesses for the FSP bounds
@@ -558,9 +560,9 @@
558560
diag(1./(log10PriorStd_a(ModelGR_a.fittingOptions.modelVarsToFit)*log(10)).^2); % Add prior in log space.
559561

560562
% for ModelGR_b parameters.
561-
ModelGR_b_fimResults = ModelGRfit_b{1}.computeFIM([],'log'); % Compute individual FIMs
563+
%ModelGR_b_fimResults = ModelGRfit_b{1}.computeFIM([],'log'); % Compute individual FIMs
562564

563-
fimGR_b_withPrior = ModelGRfit_b{1}.totalFim(ModelGR_b_fimResults,ModelGRfit_b{1}.dataSet.nCells,diag(1./(log10PriorStd_b(ModelGRfit_b{1}.fittingOptions.modelVarsToFit)*log(10)).^2)); % Add prior in log space.
565+
%fimGR_b_withPrior = ModelGRfit_b{1}.totalFim(ModelGR_b_fimResults,ModelGRfit_b{1}.dataSet.nCells,diag(1./(log10PriorStd_b(ModelGRfit_b{1}.fittingOptions.modelVarsToFit)*log(10)).^2)); % Add prior in log space.
564566

565567
%
566568
if min(eig(fimGR_a_withPrior))<1
@@ -570,12 +572,12 @@
570572
fimGR_a_covFree = fimGR_a_withPrior^-1;
571573
fimGR_a_covFree = 0.5*(fimGR_a_covFree+fimGR_a_covFree');
572574

573-
if min(eig(fimGR_b_withPrior{1}))<1
574-
disp('Warning -- FIM has one or more small eigenvalues. Reducing proposal width to 10x in those directions. MH Convergence may be slow.')
575-
fimGR_b_withPrior{1} = fimGR_b_withPrior{1} + 1*eye(length(fimGR_b_withPrior{1}));
576-
end
577-
fimGR_b_covFree = fimGR_b_withPrior{1}^-1;
578-
fimGR_b_covFree = 0.5*(fimGR_b_covFree+fimGR_b_covFree');
575+
% if min(eig(fimGR_b_withPrior{1}))<1
576+
% disp('Warning -- FIM has one or more small eigenvalues. Reducing proposal width to 10x in those directions. MH Convergence may be slow.')
577+
% fimGR_b_withPrior{1} = fimGR_b_withPrior{1} + 1*eye(length(fimGR_b_withPrior{1}));
578+
% end
579+
% fimGR_b_covFree = fimGR_b_withPrior{1}^-1;
580+
% fimGR_b_covFree = 0.5*(fimGR_b_covFree+fimGR_b_covFree');
579581

580582

581583
%% Run MH on GR Models.
@@ -693,56 +695,67 @@
693695

694696
save('conv2solnTensor_postData','GRpars_a','combinedGRModel_a', 'ModelGRfit_a','ModelGRfit_b','GRpars_b','GR_b_fspSoln', 'fspSolnsSMM')
695697

696-
dusp1 = input('(0) GR-beta has no effect on DUSP1;\n(1) GR-beta turns off the DUSP1 gene');
698+
dusp1 = input('(0) GR-beta turns off the DUSP1 gene;\n(1) GR-beta has no effect on DUSP1;\nChoose your destiny: ');
697699

698700
switch dusp1
699701
case 0
700702
%% STEP 2.A.1. -- Add DUSP1 to the existing GR model.
701703
% Copy parameters from the 100nM Dex stim case in GR.
702704
fitIters = 3;
703-
ModelGRDusp100nM = ModelGRfit_a{3};
704-
ModelGRDusp100nM.species = {'offGene';'onGene';'cytGR_a';'nucGR_a';
705-
'rna';'cytGR_b';'nucGR_b'};
706-
ModelGRDusp100nM.initialCondition = [2;0;24;1;5;12;12];
707-
ModelGRDusp100nM.propensityFunctions = {'kon*offGene*nucGR_a';'koff*onGene';
708-
'(kcn0 + (t>0)*kcn1*IDex/(MDex+IDex))*cytGR_a';'knc*nucGR_a';
709-
'kg1';'gg1*cytGR_a';'gg2*nucGR_a';'kr*onGene';'gr*rna';
710-
'kcn2*cyt_b';'knc*nucGR_b';'kg2';'gg1*cytGR_b';'gg2*nucGR_b'};
711-
ModelGRDusp100nM.stoichiometry = [-1,1,0,0,0,0,0,0,0,0,0,0,0,0;...
712-
1,-1,0,0,0,0,0,0,0,0,0,0,0,0;...
713-
0,0,-1,1,1,-1,0,0,0,0,0,0,0,0;...
714-
0,0,1,-1,0,0,-1,0,0,0,0,0,0,0;...
715-
0,0,0,0,0,0,0,1,-1,0,0,0,0,0;...
716-
0,0,0,0,0,0,0,0,0,-1,1,1,-1,0;...
717-
0,0,0,0,0,0,0,0,0,1,-1,0,0,-1];
718-
ModelGRDusp100nM.useHybrid = true;
719-
ModelGRDusp100nM.hybridOptions.upstreamODEs = {'cytGR_a','nucGR_a','cytGR_b','nucGR_b'};
720-
ModelGRDusp100nM.solutionScheme = 'FSP';
721-
ModelGRDusp100nM.customConstraintFuns = [];
722-
ModelGRDusp100nM.fspOptions.bounds = [0;0;0;2;2;2;2;400];
723-
ModelGRDusp100nM.fittingOptions.modelVarsToFit = 1:4;
724-
ModelGRDusp100nM = ModelGRDusp100nM.formPropensitiesGeneral('EricModDusp1');
725-
duspLogPrior = @(x)-sum((log10(x(:))'-log10PriorMean(1:4)).^2./(2*log10PriorStd(1:4).^2));
726-
ModelGRDusp100nM.fittingOptions.logPrior = duspLogPrior;
705+
ModelGRDusp100nM_a = ModelGRfit_a{3};
706+
ModelGRDusp100nM_b = ModelGRfit_b{1};
707+
ModelGRDusp100nM_a.species = {'cytGR_a';'nucGR_a';'offGene';'onGene';'rna'};
708+
ModelGRDusp100nM_b.species = {'cytGR_b';'nucGR_b';'offGene';'onGene';'rna'};
709+
ModelGRDusp100nM_a.initialCondition = [5;1;2;0;5];
710+
ModelGRDusp100nM_b.initialCondition = [5;1;2;0;5];
711+
ModelGRDusp100nM_a.propensityFunctions = {'(kcn0 + (t>0)*kcn1*IDex/(MDex+IDex)) * cytGR_a'; 'ba1'; 'dc * cytGR_a'; 'kn2c * nucGR_a'; 'dn * nucGR_a';...
712+
'kon*offGene*nucGR_a';'koff*onGene';'kr*onGene';'dr*rna'};
713+
ModelGRDusp100nM_b.propensityFunctions = {'kb1 * cytGR_b'; 'bb1'; 'dc * cytGR_b'; 'kn2c * nucGR_b'; 'dn * nucGR_b';...
714+
'kon*offGene';'koff*onGene*nucGR_b';'kr*onGene';'dr*rna'};
715+
ModelGRDusp100nM_a.stoichiometry = [-1,1,-1,1,0,0,0,0,0;...
716+
1,0,0,-1,-1,0,0,0,0;...
717+
0,0,0,0,0,-1,1,0,0;...
718+
0,0,0,0,0,1,-1,0,0;...
719+
0,0,0,0,0,0,0,1,-1];
720+
ModelGRDusp100nM_b.stoichiometry = [-1,1,-1,1,0,0,0,0,0;...
721+
1,0,0,-1,-1,0,0,0,0;...
722+
0,0,0,0,0,-1,1,0,0;...
723+
0,0,0,0,0,1,-1,0,0;...
724+
0,0,0,0,0,0,0,1,-1];
725+
ModelGRDusp100nM_a.useHybrid = true;
726+
ModelGRDusp100nM_b.useHybrid = true;
727+
ModelGRDusp100nM_a.hybridOptions.upstreamODEs = {'cytGR_a','nucGR_a'};
728+
ModelGRDusp100nM_b.hybridOptions.upstreamODEs = {'cytGR_b','nucGR_b'};
729+
ModelGRDusp100nM_a.solutionScheme = 'FSP';
730+
ModelGRDusp100nM_b.solutionScheme = 'FSP';
731+
ModelGRDusp100nM_a.customConstraintFuns = [];
732+
ModelGRDusp100nM_b.customConstraintFuns = [];
733+
ModelGRDusp100nM_a.fspOptions.bounds = [0;0;0;2;2;400];
734+
ModelGRDusp100nM_b.fspOptions.bounds = [0;0;0;2;2;400];
735+
ModelGRDusp100nM_a.fittingOptions.modelVarsToFit = 10:13;
736+
ModelGRDusp100nM_b.fittingOptions.modelVarsToFit = 6:9;
737+
ModelGRDusp100nM_a = ModelGRDusp100nM_a.formPropensitiesGeneral('EricModDusp1_a');
738+
ModelGRDusp100nM_b = ModelGRDusp100nM_b.formPropensitiesGeneral('EricModDusp1_b');
739+
duspLogPrior_a = @(x)-sum((log10(x(:))'-log10PriorMean(10:13)).^2./(2*log10PriorStd(10:13).^2));
740+
duspLogPrior_b = @(x)-sum((log10(x(:))'-log10PriorMean(6:9)).^2./(2*log10PriorStd(6:9).^2));
741+
ModelGRDusp100nM_a.fittingOptions.logPrior = duspLogPrior_a;
742+
ModelGRDusp100nM_b.fittingOptions.logPrior = duspLogPrior_b;
727743

728-
%% STEP 2.A.2. -- Load pre-fit parameters into model.
729-
if loadPrevious
730-
load('EricModelDusp1_MMDex','DUSP1pars')
731-
else
732744
%% Pull the DUSP1 parameters from the Model
733745
% Find the indices of the desired parameter names
734-
knc_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'knc'));
735-
kg1_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'kg1'));
736-
gg1_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'gg1'));
737-
gg2_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'gg2'));
746+
knc_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'kn2c'));
747+
kg1_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'ba1'));
748+
gg1_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'dc'));
749+
gg2_idx = find(strcmp(ModelGRDusp100nM.parameters(:,1), 'dn'));
738750

739751
% Extract the values and store them in DUSP1pars
740752
DUSP1pars = [ModelGRDusp100nM.parameters{knc_idx, 2}, ...
741753
ModelGRDusp100nM.parameters{kg1_idx, 2}, ...
742754
ModelGRDusp100nM.parameters{gg1_idx, 2}, ...
743755
ModelGRDusp100nM.parameters{gg2_idx, 2}];
744-
end
756+
745757
ModelGRDusp100nM.parameters(ModelGRDusp100nM.fittingOptions.modelVarsToFit,2) = num2cell(DUSP1pars);
758+
746759
%% STEP 2.A.3. -- Load and Associate with DUSP1 smFISH Data (100nM Dex Only)
747760
% The commented code below would be needed to fit multiple conditions,
748761
% but that is not used in this case. It is left here in case it is
@@ -876,10 +889,10 @@
876889
ModelGRDusp100nM.species = {'offGene';'onGene';'cytGR_a';'nucGR_a';
877890
'rna';'cytGR_b';'nucGR_b'};
878891
ModelGRDusp100nM.initialCondition = [2;0;24;1;5;12;12];
879-
ModelGRDusp100nM.propensityFunctions = {'kon*offGene*nucGR_a';'koff*onGene*nucGR_b';
892+
ModelGRDusp100nM.propensityFunctions = {'kon*offGene*nucGR_a';'koff*onGene';
880893
'(kcn0 + (t>0)*kcn1*IDex/(MDex+IDex))*cytGR_a';'knc*nucGR_a';
881-
'kg1';'gg1*cytGR_a';'gg2*nucGR_a';'kr*onGene';'gr*rna';
882-
'kcn2*cyt_b';'knc*nucGR_b';'kg2';'gg1*cytGR_b';'gg2*nucGR_b'};
894+
'ba1';'dc*cytGR_a';'dn*nucGR_a';'kr*onGene';'gr*rna';
895+
'kcn2*cyt_b';'kn2c*nucGR_b';'bb1';'dc*cytGR_b';'dn*nucGR_b'};
883896
ModelGRDusp100nM.stoichiometry = [-1,1,0,0,0,0,0,0,0,0,0,0,0,0;...
884897
1,-1,0,0,0,0,0,0,0,0,0,0,0,0;...
885898
0,0,-1,1,1,-1,0,0,0,0,0,0,0,0;...
@@ -892,9 +905,9 @@
892905
ModelGRDusp100nM.solutionScheme = 'FSP';
893906
ModelGRDusp100nM.customConstraintFuns = [];
894907
ModelGRDusp100nM.fspOptions.bounds = [0;0;0;2;2;2;2;400];
895-
ModelGRDusp100nM.fittingOptions.modelVarsToFit = 1:4;
908+
ModelGRDusp100nM.fittingOptions.modelVarsToFit = 10:13;
896909
ModelGRDusp100nM = ModelGRDusp100nM.formPropensitiesGeneral('EricModDusp1');
897-
duspLogPrior = @(x)-sum((log10(x(:))'-log10PriorMean(1:4)).^2./(2*log10PriorStd(1:4).^2));
910+
duspLogPrior = @(x)-sum((log10(x(:))'-log10PriorMean(10:13)).^2./(2*log10PriorStd(10:13).^2));
898911
ModelGRDusp100nM.fittingOptions.logPrior = duspLogPrior;
899912

900913
%% STEP 2.A.2. -- Load pre-fit parameters into model.

0 commit comments

Comments
 (0)