In this work we study a few basic questions for PDE learning from observed solution data. Using various types of PDEs, we show 1) how the approximate dimension (richness) of the data space spanned by all snapshots along a solution trajectory depends on the differential operator and initial data, and 2) identifiability of a differential operator from its solution data and possible instabilities. Then we propose a consistent and sparse local regression method for general PDE identification. Our method requires minimal amount of local measurements in space and time from a single solution by enforcing global consistency and sparsity.