I am still pretty new to Mathematica, but have to solve a quite complex problem.

The function

```
kdis[x_, a_] := KernelMixtureDistribution[Inner[List, x, (1+x)^a, List]]

```

gives me a bivariate kernel density estimate of two gamma random variables given by

```
gamran[k_, t_, nr_] := RandomVariate[GammaDistribution[k, t], nr]

```

Given some observations 

```
xy = {{x_1,y_1},{x_2,y_2},...,{x_n,y_n}}

```

, I seek to maximize the likelihood function

```
Likelihood[kdis[gamran[2, 0.5, 250000], gamran[k, t, 250000]], xy]

```

using either NMaximize or FindDistributionParameters. 

```
NMaximize[{Likelihood[kdis[gamran[2, 0.5, 250000], gamran[k, t, 250000]], xy], 0 < k <= 10 && 0 < t <= 10}, {{k, 0.01, 10}, {t, 0.01,10}},StepMonitor :> Print[k]]

FindDistributionParameters[xy,kdis[gamran[2.7, 0.5, 250000], gamran[k,t,250000]],{{k,1},{t,1}}, ParameterEstimator -> "MaximumLikelihood"]

```

Both functions seem to run endlessly. In particular, NMaximize never even finishes the first step (StepMonitor never prints something). As the problem is constrained to positive values of $k$ and $t$, I define constraints and a positive initial region for NMaximize. However, I get the following error message:

```
GammaDistribution::posprm: "Parameter k at position 1 in GammaDistribution[k,t] is expected to be positive."

```

which makes no sense to me given the constraints and the intial region. 

I know the problem is computationally very demanding. However, plotting the Likelihood with ContourPlot takes some time, but works and gives me relatively well defined maxima within my initial region.
So I dont quite understand what is wrong with my script.