function [llf, h] = loglike_RV_N_N(est, RV, r) % Calculates the loglikelihood assuming it follows a normal distribution
    alpha0    = est(1);                        % Sets alpha0 to first column value from estimation matrix
    alpha1    = est(2);                        % Sets alpha1 to second column value from estimation matrix
    beta1     = est(3);                        % Sets beta1 to third column value from estimation matrix
    xi        = est(4);
    phi       = est(5);
    tau1      = est(6);
    tau2      = est(7);
    sigma2u   = est(8);
    mean_ret  = est(9);
    start_vol = est(10);

    N         = length(r);                      % Set loop length N to size of return series a
    h         = zeros(N,1);
    ln_h      = zeros(N,1);                     % Fills volatility estimates vector with zeros
    epsilon   = zeros(N,1);
    ln_RV_e   = zeros(N,1);
    u         = zeros(N,1);
    a         = r - mean_ret;
    a2        = a.^2;
    ln_RV     = log(RV);
    
    % Parameter restrictions
    if alpha0 <= 0
        llf = -inf;
        return;
    end
    
    if alpha1 < 0
        llf = -inf;
        return;
    end
    
    if beta1 < 0
        llf = -inf;
        return;
    end
    
    if (phi*alpha1) + beta1 >= 1
        llf = -inf;
        return;
    end
        
    if sigma2u <= 0
        llf = -inf;
        return;
    end
    
    if start_vol > 10
        llf = -inf;
        return;
    end

    if start_vol < -10
        llf = -inf;
        return;
    end
    
    for t = 1:N
        if (t == 1)                               % If T is equal to 1. Needed as we don't have previous observations
            ln_h(t) = start_vol;                  % Sets volatility ln_h to parameter start_vol
            h(t) = exp(ln_h(t));
            epsilon(t) = a(t) ./ sqrt( h(t) );
            ln_RV_e(t) = xi + phi * ln_h(t) + tau1 * epsilon(t) + tau2 * (epsilon(t).^2 - 1);
            u(t) =  ln_RV_e(t) - ln_RV(t);
        else
            ln_h(t) = alpha0 + alpha1 * ln_RV(t-1) + beta1 * ln_h(t-1); % Calculates volatility ln_h 
            h(t) = exp(ln_h(t));
            epsilon(t) = a(t) ./ sqrt( h(t) );
            ln_RV_e(t) = xi + phi * ln_h(t) + tau1 * epsilon(t) + tau2 * (epsilon(t).^2 - 1);
            u(t) =  ln_RV_e(t) - ln_RV(t);
        end
    end 
    
    u2 = u.^2;
    llf = -0.5 .* sum(ln_h + a2./h) -0.5 .* sum(log(sigma2u) + u2./sigma2u)  + log(1/sigma2u);  
    % Loglikehood function for normal distribution for GARCH and
    % measurement equation
    % + log(1/sigma2u) Jeffery's Prior  + log(1/sigma2u)
end
    
   