function[Tran,s,probst,alambda,asigmay] = markovappr(lambda,sigma,m,N)

stvy = sqrt(sigma^2/(1-lambda^2));

ymax = m*stvy;
ymin = -ymax;

w = (ymax-ymin)/(N-1);

s = ymin:w:ymax;

for j=1:N;
    
    for k=2:N-1;
        
        Tran(j,k) = normcdf(s(k)-lambda*s(j)+w/2,0,sigma)...
            - normcdf(s(k)-lambda*s(j)-w/2,0,sigma);
    end;
    
    Tran(j,1) = normcdf(s(1)-lambda*s(j) + w/2,0,sigma);
    Tran(j,N) = 1 - normcdf(s(N) - lambda*s(j)-w/2,0,sigma);
    
end;

if sum(Tran') ~= ones(1,N);
    str = find(sum(Tran')-ones(1,N));
    disp('error in transition matrix');
    disp(['rows ',num2str(str),' does not sum to one']);
    
end;

Trans=Tran';
probst = (1/N)*ones(N,1);
test=1;

while (test > 10^(-8));
    probst1 = Trans*probst;
    test=max(abs(probst1-probst));
    probst=probst1;
end;

meanm = s*probst;
varm = ((s-meanm).^2)*probst;

midaut1 = (s-meanm)'*(s-meanm);

probmat = probst*ones(1,N);

midaut2 = Tran.*probmat.*midaut1;

autocov1=sum(sum(midaut2));

alambda=autocov1/varm;

asigmay = sqrt(varm);