Array indices must be positive integers or logical values.

1 次查看(过去 30 天)
clearvars
% Setup
t_max= 1200; % in s
k=0.01; % in s^-1
g_nr=2; % in molecules/s
t_nr=0; % in s
x_nr=5; % in molecules
time_log=0;
% Simulation of the non-regulated construct
t=0;
prod_rate=g_nr;
deg_rate=k*x_nr;
while t<t_max
wait_time=-log(rand)/(prod_rate+(deg_rate*x_nr(t_nr)));
prob_prod=prod_rate/(prod_rate+deg_rate*x_nr(t_nr)); % Propensity of production
t=t+wait_time; % Update current time
t_nr=t_nr+1; % Update the number of steps (reactions) associated with the experiment.
time_log(t_nr)=t; % Add the current time to the time log
if rand < prob_prod % Defines whether production takes place based on Monte Carlo method.
x_nr(t_nr)=x_nr(t_nr-1)+1; % Implements production
else
x_nr(t_nr)=x_nr(t_nr-1)-1; % Implements degradation
end
end
if t>t_max
t_nr(end)=[];
x_nr(end)=[];
end
close all
plot(time_log,t_nr,'r','LineWidth',2)
hold on
plot(time_log,x_nr,'r','LineWidth',2)
I get the error Array indices must be positive integers or logical values in line 24:
wait_time=-log(rand)/(prod_rate+(deg_rate*x_nr(t_nr)));
Does anyone know why this is?
  7 个评论
Jonathan Louie
Jonathan Louie 2022-6-6
I checked the units and I understand that the deg rate should be 1/s molecules/s. But the instructions given for my assignment is:
the degradation rate of both models will be given by: deg(x(t)) = kx(t) where x(t) stands for the number of repressor molecules in the system. k is given as 0.01s^-1, but what I'm unsure of is if the statement below is referring to the variable for number of repressor molecules.
Define a variable (x_nr) with the number of molecules of the non-regulated system and initialize it with value 5 molecules. As for tnr &, xnr will be a vector that contains the changes in number of molecules corresponding to the different time points of vector tnr.
Because that's why I plugged in x_nr for x(t) for the deg rate equation.
Also for reference, I'm trying to plot the dynamics of the non-regulated construct and the auto-regulated construct by running the Gillespie algorithm until we reach the maximum time, t_max. To do this, set up a while loop so that the state of the system’s dynamics – its transitions – is simulated from 0 s up to time t_max, then drawing a wait time dt from an Negative exponential distribution (Poisson process).
Thank you so much for any help I really appreciate it.
Image Analyst
Image Analyst 2022-6-6
And what about the help you got below, you know, in the official Answers section of this page? Did you see it? If not, scroll down.

请先登录,再进行评论。

回答(2 个)

Image Analyst
Image Analyst 2022-6-5
Since this is one of the most common FAQs, see the FAQ for a thorough discussion:

Torsten
Torsten 2022-6-6
if rand < prob_prod % Defines whether production takes place based on Monte Carlo method.
x_nr(t_nr)=x_nr(t_nr-1)+1; % Implements production
else
x_nr(t_nr)=x_nr(t_nr-1)-1; % Implements degradation
end
The first time you enter this if-clause, t_nr = 1.
Thus you address x_nr(t_nr-1) = x_nr(0).
This throws an error in MATLAB since array indices start with index 1, not 0.

类别

Help CenterFile Exchange 中查找有关 Scripts 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by