In this paper, we study the persistence, boundedness, convergence, invariance and global asymptotic behavior of the positive solutions of the second order difference system xn+1 = α1 + ae−xn−1 + byne−yn−1 , (0.1) yn+1 = α2 + ce−yn−1 + dzne−zn−1, zn+1 = α3 + he−zn−1 + jxne−xn−1, n = 0, 1, 2,.... Here xn, yn, zn can be considered as population densities of three species such that the population density of xn, yn, zn depends on the growth of yn, zn, xn respectively with growth rate b, d, j respectively. The positive real numbers α1, α2, α3 are immigration rate of xn, yn, zn respectively, while a, c, h denotes the growth rate of xn, yn, zn respectively, and the initial values x−1, y−1, z−1, x0, y0, z0 are nonnegative numbers.