  The current analysis aims to present the results of a QCD analysis of
diffractive parton distribution functions (diffractive PDFs) at next-to-leading
order (NLO) accuracy in perturbative QCD. In this new determination of
diffractive PDFs, we use all available and up-to-date diffractive deep
inelastic scattering (diffractive DIS) datasets from H1 and ZEUS Collaborations
at HERA including the most recent H1/ZEUS combined measurements. In this
analysis, we consider the heavy quark contributions to the diffractive DIS in
the so-called framework of {\tt FONLL} general mass variable flavor number
scheme (GM-VFNS). The uncertainties on the diffractive PDFs are calculated
using the standard "Hessian error propagation" which served to provide a more
realistic estimate of the uncertainties. This analysis are enriched, for the
first time, by including the nonperturbative higher twist (HT) effects in the
calculation of diffractive DIS cross sections which are particularly important
at large-$x$ and low $Q^{2}$ regions. Then, the stability and reliability of
the extracted diffractive PDFs are investigated upon inclusion of HT effects.
We discuss the novel aspects of the approach used in this QCD fit, namely,
optimized and flexible parameterizations of diffractive PDFs, the inclusion of
HT effects, and considering the recent H1/ZEUS combined dataset. Finally, we
present the extracted diffractive PDFs with and without the presence of HT
effects, and discuss the fit quality and the stability upon variations of the
kinematic cuts and the fitted datasets. We show that the inclusion of HT
effects in diffractive DIS can improve the description of the data which leads,
in general, to a very good agreement between data and theory predictions.
