<jupyter_start><jupyter_text><jupyter_code>import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp
%config InlineBackend.figure_format = "svg"
def growth_rate(y):
  L1 = (1-y[2])*1 + 1.5*y[2]
  L2 = (1-y[2])*1.5 + 1.*y[2]
  D = L1*y[0]/(y[0]+y[1]) + L2*y[1]/(y[0]+y[1]) 
  return D,L1,L2

def switching_ode(t, y,h):
  y_new = np.zeros(3)
  D,L1,L2 = growth_rate(y)   
  y_new[0] = (L1-D)*y[0] + h*(y[1] - y[0]) 
  y_new[1] = (L2-D)*y[1] - h*(y[1] - y[0]) 
  return y_new

def run_population(y0,environment,h):
  Y = []
  T = []
  epochs = environment[0]
  states = environment[1]
  y_last = [y0[0],y0[1],states[0]]
  
  for k in range(len(epochs))[:-1]:
    sol = solve_ivp(switching_ode, [epochs[k],epochs[k+1]], y_last,args=[h])
    Y.append(sol.y)
    T.append(sol.t)
    y_last = [sol.y[0,-1],sol.y[1,-1],states[k+1]]
  return np.concatenate(T),np.concatenate(Y,axis=1)


T,Y = run_population([0.5,1],[[0,100,200],[0,1,0]],0.004)
T,Y = T[10:-1],Y[:,10:-1]
fig,axs = plt.subplots(ncols = 2,figsize=(10,3))
 
ax = axs[0]
ax.plot(T,Y[1,:]/(Y[1,:]+Y[0,:]))
ax.plot(T,Y[0,:]/(Y[1,:]+Y[0,:]))
ax.set_xlabel("time")
ax.set_ylim(0,1)

ax = axs[1]
ax.plot(T,[growth_rate(y)[0] for y in np.transpose(Y)])
ax.set_ylim(0,2)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.2566510172143975, 'fraction_numerical': 0.07276995305164319, 'mean_word_length': 4.684444444444445, 'pattern_counts': {'":': 0, '<': 5, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '1', 'dir': 'jupyter-structured-clean-dedup', 'id': '403769', 'license': 'MIT', 'n_tokens_mistral': 701, 'n_tokens_neox': 682, 'n_words': 90, 'path': 'phenotypic_switching.ipynb', 'repo_name': 'elevien/phenotypic_variability_review', 'size': '1278'}