Vector valued autoregressive models with fractionally integrated errors are considered. The possibility of the coefficient matrix of the model having eigenvalues with absolute values equal or close to unity is included. Quadratic approximation to the log-likelihood ratios in the vicinity of auxiliary estimators of the parameters is obtained and used to make a rough identification of the approximate unit eigenvalues, including complex ones, together with their multiplicities. Using the identification thus obtained, the stationary linear combinations (cointegrating relationships) and the trends that induce the nonstationarity are identified, and Wald-type inference procedures for the parameters associated with them are constructed. As in the situation in which the errors are independent and identically distributed (i.i.d.), the limiting behaviors are nonstandard in the sense that they are neither normal nor mixed normal. In addition, the ordinary least squares procedure, which works reasonably well in the i.i.d. errors case, becomes severely handicapped to adapt itself approximately to the underlying model structure, and hence its behavior is significantly inferior in many ways to the procedures obtained here.