function dydt=f(t,x)
%   All equations and values are taken from the article:
%   "Modelling the dynamics of the yeast pheromone pathway." 
%   by Kofahl et. al


%   Rate constants, units are min^-1 or min^-1*nM^-1 depending on the order
%   of the reactions

    

    k1= 0.03;
    k2= 0.0012;
    k3= 0.6;
    k4= 0.24;
    k5= 0.024;
    k6= 0.0036;
    k7= 0.24;
    k8= 0.33; %0.33 kofahl
    k9= 2000;
    k10= 0.1;
    k11= 5;
    k12= 1;
    k13= 3;
    k14= 1;
    k15= 3;
    k16= 3;
    k17= 100;
    k18= 5;
    k19= 1;
    k20= 10;
    k21= 5;
    k22= 47;
    k23= 5;
    k24= 345;
    k25= 5;
    k26= 50;
    k27= 5;
    k28= 140;
    k29= 10;
    k30= 1;
    k31= 250;
    k32= 5;
    k33= 50;
    k34= 18;
    k35= 10;
    k36= 0.1;
    k37= 0.1;
    k38= 0.01;
    k39= 18;
    k40= 1;
    k41= 0.002; %0.002 kofahl
    k42= 0.1;
    k43= 0.01;
    k44= 0.01;
    k45= 0.1;
    k46= 200;
    k47= 1;
    k48= 7.77; %Msg2 (Dyjack)
    k49= 1; %Add TEFp production rate!
    k50= 0.000005;
    k51= 0.01;
    k52= 0.1;

%   Reaction kinetics for each reaction
    v1 = x(1)*x(29)*k1;     %1
    v2 = x(2)*x(1)*k2;      %2
    v3 = x(3)*k3;           %3
    v4 = x(3)*k4;           %4
    v5 = x(2)*k5;           %5
    v6 = x(3)*x(4)*k6;      %6
    v7 = x(5)*k7;           %7
    v8 = x(5)*x(36)*k8;      %8
    v9 = x(7)*x(6)*k9;      %9
    v10 = x(6)*x(8)*k10;    %10
    v11 = x(9)*k11;         %11
    v12 = x(10)*x(11)*k12;    %12
    v13 = x(12)*k13;         %13
    v14 = x(13)*x(14)*k14;   %14
    v15 = x(15)*k15;         %15
    v16 = x(12)*x(15)*k16;   %16
    v17 = x(18)*k17;         %17
    v18 = x(9)*x(16)*k18;   %18
    v19 = x(17)*k19;         %19
    v20 = x(17)*k20;         %20
    v21 = x(17)*k21;         %21
    v22 = x(18)*k22;         %22
    v23 = x(18)*k23;         %23
    v24 = x(19)*k24;         %24
    v25 = x(19)*k25;         %25
    v26 = x(20)*k26;         %26
    v27 = x(20)*k27;         %27
    v28 = x(21)*k28;         %28
    v29 = x(22)*x(14)*k29;   %29
    v30 = x(24)*k30;         %30
    v31 = x(24)*k31;         %31
    v32 = x(22)*k32;         %32
    v33 = x(23)*k33;         %33
    v34 = x(25)*x(23)*k34;   %34
    v35 = x(26)*k35;         %35
    v36 = x(26)*x(27)*k36;   %36
    v37 = x(28)*k37;         %37
    v38 = x(28)*k38;         %38
    v39 = x(30)*((x(23)^(2))/((100^(2))+x(23)^(2)))*k39; %39
    v40 = x(31)*k40;         %40
    v41 = x(30)*x(35)*k41;   %41
    v42 = x(6)*x(31)*k42;    %42
    v43 = x(33)*k43;         %43
    v44 = x(34)*k44;         %44
    v45 = x(31)*x(35)*k45;   %45
    v46 = ((x(23)^(2))/((4^(2))+x(23)^(2)))*k46; %46
    v47 = x(36)*k47;          %47
    v48 = ((x(26)^2.1))/((100^(2.1))+x(26)^(2.1))*k48; %48
    
    
    exposure=x(37)*t;
    if exposure>5000
    v49 = k49;
    else
    v49 = 0;
    end
    
    v50 = k50*x(37)*x(40);
    v51 = k51*x(38)^2*x(39)/500;
    v52 = k52*x(41);

%   Rate equations for each specie
dydt=[-v1     %alpha factor 
      -v2+v3-v5   %Ste2 
      v2-v3-v4    %Ste2active
      -v6+v9        %Galphabetagamma
      v6-v7-v8    %GalphaGTP 
      v6-v9-v10+v11+v21+v23+v25+v27+v32-v42+v43   %Gbetagamma
      v7+v8-v9    %GaplhaGDP
      -v10+v11+v16-v17  %C 
      v10-v11-v18+v19   %D   
      -v12+v13+v17+v21+v23+v25+v27+v32  %Ste5  
      -v12+v13+v17+v21+v23+v25+v27+v32  %Ste11
      v12-v13-v16 %A
      -v14+v15+v17+v21+v23+v25+v27+v32  %Ste7 
      -v14+v15+v17+v21+v23+v25+v27-v29+v30+v33  %Fus3
      v14-v15-v16 %B
      -v18+v19+v21+v23+v25+v27+v32    %Ste20    
      v18-v19-v20-v21 %E 
      v20-v22-v23 %F 
      v22-v24-v25 %G 
      v24-v26-v27 %H 
      v26-v28+v31 %I 
      v28-v29+v30-v32   %L
      v28-v33-v34+v35   %Fus3PP
      v29-v30-v31 %K  
      -v34+v35  %Ste12 
      v34-v35   %Ste12active 
      -v36+v37  %Bar1 
      v36-v37-v38 %Bar1active
      v38 %Bar1activeex
      -v39+v40-v41    %Far1 
      v39-v40-v42+v43+v44-v45   %Far1PP
      v41 %Far1ubi
      v42-v43   %M
      -v44+v45  %N
      v44-v45   %Cdc28 
      v46-v47       %Sst2active
      v48 %Cre/Fus1
      v49 %Cas9
      0.5*v50 %Mated diploids
      -v50 %Haploids
      v51-v52 %Cas9-gRNA complex
      ];

end