A recently proposed method implementing a pseudo absorption term in the neutron transport equation for mode eigenvalue calculations in the Monte Carlo technique has been discussed. This method is known to overcome the difficulty in mode eigenvalue calculations for fissionable systems with large subcriticality. This paper has demonstrated that this technique certainly can provide stable mode eigenvalue calculations for problems that the original MCNP 4C fails to solve and that the figure of merit of the calculated value is largely improved. However, caution should be taken that the variance of calculated values is underestimated due to the inter-cycle correlation of the values. A method that provides eigenfunctions with higher order criticality eigenvalues has been applied to the second order mode eigenvalue calculations. The method partitions a whole space into two regions. The estimate of the multiplication factor in each region is forced to be equal to each other. This method is found to be applicable to calculations for the second mode eigenfunction. The conditions on convergence of the second mode eigenfunction are discussed. Hotelling's method, which explicitly subtracts lower order eigenfunctions from the fission source distribution, has been applied to a test problem. The third and fourth order forward eigenfunctions are successfully obtained with the use of Hotelling's method.