Incorporating a natural death component into a Monod equation

Question

Hi I have a series of monod equations representing a system with two substrates and two species. Each of the species is producing a substrate that the other species consumes as such the success of species one is dependent on species two and vice versa. I can simulate these results in python through a series of ODEs however when I attempt to incorporate a death rate into the species I am having a hard time getting results which are not negative (biologically impossible).

I've read that I should be incorporating a death rate constant and then just multiplying by the number of a species present like so

Sp1=GM1*(y[3]/(ks1+y[3])*y[0])- kd*y[0]  

Where; GM1 = max growth rate of sp1 ks1 = concentration of substrate when usage is at 1/2 GM1 kd=cell death rate y[3] is the substrate being utilised y[0] is the species

however this approach is causing the species numbers to become negative.

I've attached the functioning python code without the cell death for anyone who would like to take a stab at it/examine the system of equations

GM1=5*10**-2 #Growth Rate Sp1    
ks1=5*10**1  #Concentration of substrate 1 at rate 1/2 GM1 max 
GM2=4*10**-3  #Growth Rate Sp2
ks2=4*10**1  #Concentration of substrate 2 at rate 1/2 GM2 max


# Conversion of substrates into biomass and creation of substrates by biomass

con=40    
con2=90      
con3=7       
con4=60  




def chemanalysis (y,t):
    Sp1=GM*(y[3]/(ks1+y[3])*y[0])  

    Sp2=GM2*(y[2]/(ks2+y[2])*y[1])  

    Sub1=1/con3*GM*(y[3]/(ks1+y[3])*y[0])  - (1/con4)* (y[2]/(ks2+y[2])*y[1])

    Sub2=1/con*GM3*(y[2]/(ks2+y[2])*y[1])  - (1/con2)* (y[3]/(ks1+y[3])*y[0])


    return [Sp1,Sp2,Sub1,Sub2]

    y0=(20,10,50,30,20,20,20,20)
    tspan=np.arange(0,10000,1)

    Conc= odeint(chemanalysis,y0,tspan )

    ALGC                = Conc[:,0]
    METC                = Conc[:,1]
    Cin                 = Conc[:,2]
    Cout                = Conc[:,3]


plt.plot(tspan,ALGC,label='Algae')
plt.plot(tspan,METC,label='Methanogens')
plt.plot(tspan,Cin,label='C in')
plt.plot(tspan,Cout,label='C out')

plt.legend()
p.show()

Any help would be greatly appreciated

Thanks

Answer 1

Do what you've read: incorporate a death rate constant and then just multiply the number of a species present by it.

In more detail:

1. Calculate the numbers the way you were doing it before, without any death rate;

   Sp1=GM1*(y[3]/(ks1+y[3])*y[0])

2. Now multiply Sp1 by 0.99 or something.

   Sp1=0.99*Sp1

This way you'll never run into negative numbers.