nonlinear fitting experimental data

Good afternoon, sorry, I'm new to matlab and have many doubts, I tried searching previous threads if you have talked about something, and so it is, but I have not found any solution that works me ... my problem is that I need to adjust some experimental two differential equations data, I worked on a code that collects data from one of them but I can not fit the other, I get error all the time ... I leave the code in case they can serve as a guide ... and not try to do that ... but I can not move ... thanks in advance.

cb exp and F2 appearing in f_ej3, is what I want to add, but do not know how, I do not really know which method to use, I have read fminsearch, RK4...but i don't know ... Can you help me?

Prueba.m

% Fichero principal. Optimizacion
% Definimos un vector inicial, pero teniendo en cuenta que este vector inicial
% es solo un conjunto inicial de valores de K0,E y N para probar
X0 = [10 2000 3 5 3000 1 2];
X = fminsearch('f_ej2',X0);
% Con los valores obtenidos representamos para ver la bondad del ajuste
% Introducimos los valores experimentales
Tiempo = [0 0.5 1 1.5 2 2.5 3 3.5 4 5 6 7 8];
CAexp = [0.0406 0.0331 0.0268 0.0208 0.0151 0.0106 0.0071 0.0057 0.0045 0.0000 0.0000 0.0000 0.0000]';
%Cbexp = [0.0607 0.0794 0.1005 0.1115 0.1150 0.1180 0.1266 0.1278 0.1268 0.1262 0.1258 0.1242 0.1184]';

% Reasignamos las constantes, pero antes las definimos como globales, para que puedan ser
% utilizadas por el fichero funcion
global K0 E N K2 E2 N2 K1
K0 = X(1);
E = X(2);
N = X(3);
K2 = X(4);
E2 = X(5);
N2 = X(6);
K1 = X(7);

% Preparando la aplicacion del ODE113
PuntoInicial = 0;
PuntoFinal = 25;
Ca0 = 0.26852;

[T,CAteo] = ode113('f_ej3',Tiempo,Ca0);
%[T1,CBteo] = ode113('f_ej3',Tiempo,Cb0);
% Hacemos la representacion
figure(1)
plot(Tiempo,CAexp,Tiempo,CAteo)
title('Ejemplo de optimizacion')
xlabel('Tiempo (s)')
ylabel('Ca (mol/L)')
legend('CAexp','CAteo')
grid

function E=f_ej2(X)
% Fichero de funcion del ejemplo 2. Minimizacion
% X(1) Corresponde a la constante K0
% X(2) Corresponde a la constante E
% X(3) Corresponde a la constante n
% Introducimos los valores experimentales
Tiempo = [0 0.5 1 1.5 2 2.5 3 3.5 4 5 6 7 8];
CAexp = [0.0406 0.0331 0.0268 0.0208 0.0151 0.0106 0.0071 0.0057 0.0045 0.0000 0.0000 0.0000 0.0000]';
%Cbexp = [0.0607 0.0794 0.1005 0.1115 0.1150 0.1180 0.1266 0.1278 0.1268 0.1262 0.1258 0.1242 0.1184]';

% Reasignamos las constantes, pero antes las definimos como globales, para que puedan ser
% utilizadas por el fichero funcion
global K0 E1 N K2 E2 N2 K1
K0 = X(1);
E1 = X(2);
N = X(3);
K2 = X(4);
E2 = X(5);
N2 = X(6);
K1 = X(7);

% Preparando la aplicacion del ODE113
PuntoInicial = 0;
PuntoFinal = 25;
Ca0=0.26852;
%Cb0 = 0;
[T,CAteo]=ode113('f_ej3',Tiempo,Ca0);
%[T1,CBteo] = ode113('f_ej3',Tiempo,Cb0);
%Calculamos un vector error como diferencia de todos los valores al cuadrado
E = sum((Caexp-CAteo).^2);
%Error = sum((Cbexp-CBteo).^2);

end

function F = f_ej3(t,Ca)
% Fichero de funcion del ejemplo 2. Integracion
% Definimos las variables como globales, para que puedan ser utilizadas por el fichero funcion

    global K0 E N K2 E2 N2 K1
    F = (-K0*exp(-E/((273+140)*8.314))*Ca.^N - K2*exp(-E2/((273+140)*8.314))*(Ca.^N2)*(Cb.^N3)/((1+K1*Ca)^1); % i change, put (Cb^N3) but my first program it is not here.
  %F2=(K0*exp(-E/((273+140)*8.314))*Ca.^N - K2*exp(-E2/((273+140)*8.314))*(Ca.^N2)(Cb.^N3))/((1+K1*Ca)^1)-K3*(Cb.^N4);

    end

thanks