This work investigates upon the determinants of students’ literacy in mathematics and reading in Italy by analyzing longitudinal data collected in four waves 2003, 2006, 2009 and 2012 by the Program for International Student Assessment (PISA) ran by OECD. PISA data provides five plausible values in reading and and five in mathematics for each student, thus three levels of analysis are detectable in the data structure: at level-1 ten students’ plausible values (five for math and five for reading) ; at level-2 students and at level-3 schools. The main aim of our analysis is: (i) to assess the relationships between students' skills in mathematics and reading and socio-cultural background; (ii) to suggest value-added measures of students’ performances at regional level in the time span considered in the analysis; (iii) to detect trajectories of performances and to advance a system of indicators at regional level. A three Level Multilevel Regression Model with Bivariate Latent Structure at level-2 (student-level) and heteroscedastic random terms at level-3 (school-level) has been adopted to single out the factors which seem to account for the greatest variability in students' skills as well as to jointly model the unobserved heterogeneity across schools and regions.