Spontaneous breaking of Lorentz symmetry at energies on the order of the Planck energy or lower is predicted by many quantum gravity theories, implying non-trivial dispersion relations for the photon in vacuum. Consequently, gamma-rays of different energies, emitted simultaneously from astrophysical sources, could accumulate measurable differences in their time of flight until they reach the Earth. Such tests have been carried out in the past using fast variations of gamma-ray flux from pulsars, and more recently from active galactic nuclei and gamma-ray bursts. We present new constraints studying the gamma-ray emission of the galactic Crab Pulsar, recently observed up to TeV energies by the Major Atmospheric Gamma-ray Imaging Cherenkov (MAGIC) collaboration. A profile likelihood analysis of pulsar events reconstructed for energies above 400 GeV finds no significant variation in arrival time as their energy increases. Ninety-five percent CL limits are obtained on the effective Lorentz invariance violating energy scale at the level of {E}{{QG}1}> 5.5× {10}17 {GeV} (4.5× {10}17 {GeV}) for a linear, and {E}{{QG}2}> 5.9× {10}10 {GeV} (5.3× {10}10 {GeV}) for a quadratic scenario, for the subluminal and the superluminal cases, respectively. A substantial part of this study is dedicated to calibration of the test statistic, with respect to bias and coverage properties. Moreover, the limits take into account systematic uncertainties, which are found to worsen the statistical limits by about 36%-42%. Our constraints would have been much more stringent if the intrinsic pulse shape of the pulsar between 200 GeV and 400 GeV was understood in sufficient detail and allowed inclusion of events well below 400 GeV.
