Fit experimental data to find unknown parameter in system of ODE

Question

I have a system of ode including 3 diff equation. I found two parameters Tcell and CD19 by experiment over time. Here is my dataset:

data1 = {{1, 1, 6734}, {1, 1.5, 8734}, {1, 2, 12826}, {1, 2.5, 
20826}, {1, 3, 34006}, {1, 3.5, 30006}, {1, 4, 26975.5}, {1, 4.5, 
20975.5}, {1, 5, 15307}, {1, 5.5, 13307}, {1, 6, 12307}, {2, 1, 
9222.5}, {2, 1.5, 9222.5}, {2, 2, 8653}, {2, 2.5, 7523}, {2, 3, 
6002}, {2, 3.5, 5002}, {2, 4, 3996}, {2, 4.5, 3002}, {2, 5, 
1985}, {2, 6, 1237}};

now I have some known rate parameters such as kpCD19 & kaCD19 and by parametricNDSolveValue, I am trying to find unknown parameters such as kD, kkill, kdecayT, koff, kon. I also have a set of differential equations that I want to fit data to these equatsions

kpCD19 = 0.5; kaCD19 = 0.2;tmax = 6;
ode = {Tcell'[t] == -kon * Tcell[t] *CD19[t] + koff * complex[t] + 
 kD* Tcell[t] -  kdecayT Tcell[t], CD19'[t] == 
kpCD19* CD19[t] - kaCD19* CD19[t] - kon * Tcell[t] *CD19 [t] + 
 koff *complex[t] - kkill*CD19[t], complex'[t] == 
kon *Tcell[t]*CD19[t] - koff *complex[t] - kD* Tcell[t], Tcell[0] == 5000, CD19[0] == 10000, complex[0] == 0};

paramSOL = ParametricNDSolveValue[ode, {Tcell, CD19, complex}{t, 0, tmax}, {kD, kkill, kdecayT,nkoff, kon}, WorkingPrecision ->100]



model[kD_, kkill_, kdecayT_, koff_, kon_][i_, t_] := 
Through[paramSOL[kD, kkill, kdecayT, koff, kon][t], List][[i]] /; And @@ NumericQ /@ {kD, kkill, kdecayT, koff, kon, i, t}    

then I fit my model by NonlinearModelfit but it seems that it did not fit it good since the Rsquared is 0.09

fitted = NonlinearModelFit[data1, model[kD, kkill, kdecayT, koff, kon][i, t], {{kD, 0.1}, {kkill, 0.1}, {kdecayT, 0.1}, {koff,   0.1}, , {kon, 0.1}}, {i, t}, WorkingPrecision -> 100] // Quiet;

fitted["RSquared"]]
fitted["ParameterTable"]
fitted["ParameterErrors"][[2]]

dataTcell = Take[data1, 12][[All, 2 ;; 3]]
dataBcell = Drop[data1, 12][[All, 2 ;; 3]]
Show[ListPlot[{dataTcell, dataBcell}, PlotLegends -> {"Tcell"},Frame -> True], Plot[{fitted[1, t], fitted[2, t]}, {t, 0, tmax}]]

when I graph the fit also it seems that it could not fit the data correctly. My question is what is the problem that it could not fit data correctly?