The pT-differential non-linear flow modes, v4,22, v5,32, v6,33 and v6,222 for π±, K±, KS0 , p + p ¯ , Λ + Λ ¯ and ϕ-meson have been measured for the first time at sNN = 5.02 TeV in Pb-Pb collisions with the ALICE detector at the Large Hadron Collider. The results were obtained with a multi-particle technique, correlating the identified hadrons with reference charged particles from a different pseudorapidity region. These non-linear observables probe the contribution from the second and third order initial spatial anisotropy coefficients to higher flow harmonics. All the characteristic features observed in previous pT-differential anisotropic flow measurements for various particle species are also present in the non-linear flow modes, i.e. increase of magnitude with increasing centrality percentile, mass ordering at low pT and particle type grouping in the intermediate pT range. Hydrodynamical calculations (iEBE-VISHNU) that use different initial conditions and values of shear and bulk viscosity to entropy density ratios are confronted with the data at low transverse momenta. These calculations exhibit a better agreement with the anisotropic flow coefficients than the non-linear flow modes. These observations indicate that non-linear flow modes can provide additional discriminatory power in the study of initial conditions as well as new stringent constraints to hydrodynamical calculations. [Figure not available: see fulltext.].