We investigate the relative importance of initial and final state effects on azimuthal correlations of gluons in low and high multiplicity p+Pb collisions. To achieve this, we couple Yang-Mills dynamics of pre-equilibrium gluon fields (IP-GLASMA) to a perturbative QCD based parton cascade for the final state evolution (BAMPS) on an event-by-event basis. We find that signatures of both the initial state correlations and final state interactions are seen in azimuthal correlation observables, such as v2 {2P C} (pT ), their strength depending on the event multiplicity and transverse momentum. Initial state correlations dominate v2 {2P C} (pT ) in low multiplicity events for transverse momenta pT > 2 GeV. While final state interactions are dominant in high multiplicity events, initial state correlations affect v2 {2P C} (pT ) for pT > 2 GeV as well as the pT integrated v2 {2P C}.Introduction. The measured azimuthal momentum anisotropies of produced particles in heavy ion collisions are well described in the framework of event-by-event hydrodynamics. In this picture a fluctuating initial geometry, dominated by fluctuating nucleon positions in the incoming nuclei, is converted into anisotropic momentum space distributions by the pressure driven final state evolution. Hydrodynamic simulations agree well with a wide range of experimental observables from the Relativistic Heavy Ion Collider (RHIC) at Brookhaven National Laboratory and the Large Hadron Collider (LHC) at CERN [1][2][3][4].Measurements in smaller collision systems such as p+p and p+A [5], in particular those of anisotropies in multiparticle correlation functions, have shown very similar features as those in heavy ion collisions. While calculations within the hydrodynamic framework have been quite successful in describing observables in these small collision systems, alternative explanations relying entirely on intrinsic momentum correlations of the produced particles can also reproduce many features of the experimental data. This includes two and more particle azimuthal correlations and their p T dependence [5][6][7] and mass splitting of identified particle v n [8]. Apart from the existence of alternative explanations, the applicability of hydrodynamics becomes increasingly doubtful as the system size decreases and gradients increase. Some recent studies argue that hydrodynamics should be applicable in systems of sizes down to ∼ 0.15 fm [9], but off-equilibrium corrections to particle distribution functions for momenta p T 0.5 GeV can be significant [10], which limits at least the quantitative reliability of the framework.So far all calculations of multi-particle correlations in small collision systems have studied either only intrinsic momentum correlations or purely final state driven ef-