Context : This work is part of a systematic X-ray survey of the Taurus star forming complex with XMM-Newton .

Aims : We study the time series of all X-ray sources associated with Taurus members , to statistically characterize their X-ray variability , and compare the results to those for pre-main sequence stars in the Orion Nebula Cluster and to expectations arising from a model where all the X-ray emission is the result of a large number of stochastically occurring flares .

Methods : The analysis of the lightcurves is based on a maximum likelihood algorithm that segments the time series in intervals of constant signal without the need of binning .
Flares are defined with criteria that take into account the amplitude and the derivative of the segmented lightcurves .
Variability statistics are evaluated for different classes of pre-main sequence stars ( protostars , cTTS , wTTS , brown dwarfs ) , and for different spectral type ranges .
Flare frequency and energy distribution are computed .

Results : We find that roughly half of the detected X-ray sources show variability above our sensitivity limit , and in \sim 26 % of the cases this variability is recognized as flares .
Variability is more frequently detected at hard than at soft energies .
The variability statistics of cTTS and wTTS are undistinguishable , suggesting a common ( coronal ) origin for their X-ray emission .
The frequency of large flares ( E > 10 ^ { 35 } erg ) on Taurus members is 1 event per star in 800 ksec .
The typical duration of these flares – probably biased by the finite observing time – is about 10 ksec .
We have for the first time applied a rigorous maximum likelihood method in the analysis of the number distribution of flare energies on pre-main sequence stars .
In its differential form this distribution follows a power-law with index \alpha = 2.4 \pm 0.5 , in the range typically observed on late-type stars and the Sun .

Conclusions : The signature of the X-ray variability in the pre-main sequence stars in Taurus and Orion provides twofold support for coronal heating by flares : ( i ) The correlation between the maximum variability amplitude and the minimum emission level indicates that both flare and quiescent emission are closely related to the coronal heating process .
( ii ) The power-law index \alpha derived for the flare energy distribution is large enough to explain the heating of stellar coronae by nano-flares ( \alpha > 2 ) , albeit associated with a rather large uncertainty that leaves some doubt on this conclusion .