In this paper, we develop an efficient lattice Boltzmann (LB) model for simulating immiscible incompressible $N$-phase flows $(N geq 2)$ based on the Cahn-Hilliard phase field theory. In order to facilitate the design of LB model and reduce the calculation of the gradient term, the governing equations of the $N$-phase system are reformulated, and they satisfy the conservation of mass, momentum and the second law of thermodynamics. In the present model, $(N-1)$ LB equations are employed to capture the interface, and another LB equation is used to solve the Navier-Stokes (N-S) equations, where a new distribution function for the total force is delicately designed to reduce the calculation of the gradient term. The developed model is first validated by two classical benchmark problems, including the tests of static droplets and the spreading of a liquid lens, the simulation results show that the current LB model is accurate and efficient for simulating incompressible $N$-phase fluid flows. To further demonstrate the capability of the LB model, two numerical simulations, including dynamics of droplet collision for four fluid phases and dynamics of droplets and interfaces for five fluid phases, are performed to test the developed model. The results show that the present model can successfully handle complex interactions among $N$ ($N geq 2$) immiscible incompressible flows.