[R] Optimization problem with constraints

[Matthias Graser]

I'm trying to da an optimization for the followig function

Zwischenwert <- function (x)
{
	lambda<-x[1];
	mu<-x[2];
	gammal<-x[3];
	mud<-x[4];
	gammad<-x[5];
	Mittelwert <-0;
	for(t in 0:(T-1))
	{
		for(i in 0:(n-1))
		{
			for(j in i:(n-1))
			{
				Mittelwert <- Mittelwert +( phi[i+1,j+1,t+1]*((j-i)*log(lambda)-log(factorial(j-i))-lambda));
			};
			if(t>0){Mittelwert <- Mittelwert +(phisum[i+1,t+1]*(-0.5*log(2*pi)-0.5*log(t*gammal+i*gammad)-((Y[t+1]-i*mud-t*mu-0.5*t*gammal+lambda*exp(mud+0.5*gammad)-t*lambda)^2 )/(2*(t*gammal+i*gammad))));};
		};
	};
};


and it's gradient

ZWGrad <- function(x)
{
	lambda<-x[1];
	mu<-x[2];
	gammal<-x[3];
	mud<-x[4];
	gammad<-x[5];
	dlambda<-0;
	dmu<-0;
	dgamma<-0;
	dmud<-0;
	dgammad<-0;
	for (t in 0:(T-1))
	{
		for (i in 0:(n-1))
		{
			for (j in i:(n-1))
			{
				dlambda<-dlambda-phi[i+1,j+1,t+1]*(-1+(j-i)/lambda);
			};
		};
	};
	for (t in 1:(T-1))
	{
		for (i in 0:(n-1))
		{
			dlambda<-dlambda+phisum[i+1,t+1]*(t-t*lambda*exp(mud+0.5*gammad))*(Y[t+1]-i*mud-mu*t-0.5*gammal*t+lambda*t*exp(mud+0.5*gammad)-lambda*t)/(2*(gammal*t+gammad*i));
			dmu<-dmu+phisum[i+1,t+1]*t*(Y[t+1]-i*mud-t*mu-0.5*gammal*t+lambda*t*exp(mud+0.5*gammad)-lambda*t)/(gammal*t+gammad*i);
			dgamma<-dgamma+phisum[i+1,t+1]*( t*(0.5*i*gammad+Y[t+1]-i*mud-t*mu+lambda*t*exp(mud+0.5*gammad)-lambda*t)^2/(2*(gammal*t+gammad*i)^2)-t/(2*(gammal*t+gammad*i))-0.125*t);
			dmud<-dmud+phisum[i+1,t+1]*(i-lambda*t*exp(mud+0.5*gammad))*(Y[t+1]-i*mud-t*mu-0.5*gammal*t+lambda*t*exp(mud+0.5*gammad)-lambda*t)/(gammal*t+gammad*i);
			dgammad<-dgammad+phisum[i+1,t+1]*(Y[t+1]-i*mud-t*mu-0.5*gammal*t+lambda*t*exp(mud+0.5*gammad)-lambda*t)*(i*(Y[t+1]-i*mud-t*mu-0.5*gammal*t-lambda*t)-lambda*t*exp(mud+0.5*gammad)*(gammal*t+gammad*i-i))/(2*(gammal*t+gammad*i)^2);
		};
	};
	GradZW <- c(dlambda,dmu,dgamma,dmud,dgammad);
};


lambda, gammal und gammd need to be larger than 0, so I've tried different constrained methods like L-BFGS-B via optim and constrOptim, but both seem to have problems with my functions. The standard optim-method works, but sometimes it returns negative values for lambda, gammal and gammad. 
Those phi and phisums are known at runtime, n and T can be pretty large, so maybe R simple can't handel it.


Many thanks
Matthias


--