A method of estimating the MFB50 combustion index, namely the crank angle at which 50% of the fuel mass has been burnt in a cylinder of an internal combustion engine provided with a drive shaft coupled to at least one pair of position sensors (5,5*; 5A, 5B), each of which being arranged at a respective end of the drive shaft; the estimation method comprising the steps of: capturing the signals coming from the two position sensors (5,5*; 5A, 5B); determining the angular torsion (Δθ) of the drive shaft based on the signals coming from the two position sensors (5,5*; 5A, 5B); determining, by means of a frequency analysis of the angular torsion (Δθ) of the drive shaft, a first harmonic (Ω 1) of the angular torsion (Δθ) signal characterized by its own module (|Ω 1 |) and by its own phase (Arg(Ω 1)); and estimating, during normal operation,the MFB50 combustion index of each cylinder based on the phase (Arg(Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal and a first correlation coefficient (P i) of each cylinder between the MFB50 combustion index of each cylinder and the phase (Arg(Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; wherein the first correlation coefficient (P i) is determined for each cylinder, during a preliminary set-up phase, as follows: P i = MFB 50 i / Arg Ω 1 i with i = 1... n MFB50 i MFB50 combustion index based on the cylinder; Arg (Ω 1) 1 phase Arg (Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; P i first correlation coefficient between the MFB50 combustion index and the phase Arg (Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; and n number of cylinders. capturing the signals coming from the two position sensors (5,5*; 5A, 5B); determining the angular torsion (Δθ) of the drive shaft based on the signals coming from the two position sensors (5,5*; 5A, 5B); determining, by means of a frequency analysis of the angular torsion (Δθ) of the drive shaft, a first harmonic (Ω 1) of the angular torsion (Δθ) signal characterized by its own module (|Ω 1 |) and by its own phase (Arg(Ω 1)); and estimating, during normal operation,the MFB50 combustion index of each cylinder based on the phase (Arg(Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal and a first correlation coefficient (P i) of each cylinder between the MFB50 combustion index of each cylinder and the phase (Arg(Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; wherein the first correlation coefficient (P i) is determined for each cylinder, during a preliminary set-up phase, as follows: P i = MFB 50 i / Arg Ω 1 i with i = 1... n MFB50 i MFB50 combustion index based on the cylinder; Arg (Ω 1) 1 phase Arg (Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; P i first correlation coefficient between the MFB50 combustion index and the phase Arg (Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; and n number of cylinders. MFB50 i MFB50 combustion index based on the cylinder; Arg (Ω 1) 1 phase Arg (Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; P i first correlation coefficient between the MFB50 combustion index and the phase Arg (Ω 1) of the first harmonic (Ω 1) of the angular torsion (Δθ) signal; and n number of cylinders.