Performance in parameter estimation from ParametricNDSolve using varied initial conditions

I have 250 data points from a timecourse exeriment in a list, with columns specifying (1) time, (2-4) initial conditions, (5) absorbance reading. I want to fit 4 parameters (k1, k2, k3, k4) in a DAE system of 7 species to this dataset, where species7 is absorbance.

At the moment, the program needs 20 seconds for two NMinimize iterations and seems to break down after 10 minutes. I'm also using a different program (Copasi), which is able to solve the problem and perform the computations about 2 magnitudes faster. So where does this code get stuck and how do I solve it? I suspect part of the problem is that the system must be calculated for every single time point.

func[ic1,ic2,ic3] :=
species7[k1, k2, k3, k4] /. ParametricNDSolve[
{"species1-7 DAEs, species1-7 initial conditions, partly given by ic1-3"},
species7,
{t, 0, 120},
{k1, k2, k3, k4}
]

fit = NMinimize[
{
Sum
[
(data[[i, 5]] -
func[data[[i, 2]], data[[i, 3]], data[[i, 4]]][data[[i, 1]]]
)^2, {i, Length[data]}
],
{{10^-2 < k1 < 10^2},
{10^-2 < k2 < 10^2},
{10^-2 < K3 < 10^2},
{10^-2 < k4 < 10^2}}
},
{k1, k2, k3, k4},
]


Thanks

edit: here is a notebook with data https://www.dropbox.com/s/4mvnn4m8ldy7373/2013-05-16%20Stackexchange.nb please note, I'm not expecting a great fit at the moment. It takes ~1min to run

-
More people might look at this if you provide a full example and some data to work with. –  Szabolcs May 15 '13 at 22:56
This is impossible to address without the full code. –  user21 May 16 '13 at 6:53
I'm not sure I understand the problem. It takes 61 sec. to compute fit on my i7 Mac. If you solve the diff. eq. with ic1, ic2, ic3 as parameters, computing fit takes about 10 sec. (With[{sol = ParametricNDSolve[..., {k1, k2, k3, k4, ic1, ic2, ic3}]}, func[ic1_, ic2_, ic3_] := g[k1, k2, k3, k4, ic1, ic2, ic3] /. sol], assuming g in the linked notebook is species7). –  Michael E2 May 17 '13 at 1:20
I don't know why, but introducing With helped this program a lot, thanks! I'm now trying to use NonlinearModelFit instead of NMinimize to get useful statistics. Unfortunately these lines are failing (even though it seems to start): fit = Monitor[NonlinearModelFit[ data, func[m, n, o][p], {k1,k2,k3,k4}, {m, n, o, p}, Method -> NMinimize], {k1,k2,k3,k4}] // Timing NonlinearModelFit::nrnum: The function value [...] is not a real number at [...] Is it possible, in principle, to use NonlinearModelFit with this interpolating function? –  ther May 19 '13 at 13:45