This work proposes a coupled phase-field framework for the simulation of self-healing materials, integrating a probabilistic description of the chemical healing kinetics. The mechanical behavior is formulated using a degraded elasticity tensor dependent on a damage function D(s, h). This damage function is governed by two interacting thermodynamic fields: an entropy-like field s, representing the irreversible degradation of the material and crack topology, and an enthalpy-like field h, representing the availability of chemical reactants required for the healing phase. The novelty of the approach lies in the constitutive coupling of these fields. While standard phase-field damage models restrict the damage variable to a normalized range, the proposed formulation extends the upper bound of the entropy field to 1+h. Consequently, the evolution of s is driven by mechanical energy release, while the evolution of h is driven by the entropy state s (Model 1) or the elasticity energy (Model 2), mimicking the activation of self-healing agents upon damage initiation. This bidirectional coupling allows for the autonomous restoration of local stiffness properties. The resulting system of coupled non-linear equations is solved using a staggered finite element scheme. This algorithm sequentially resolves the mechanical equilibrium, the entropy evolution, and the enthalpy reaction-diffusion equation at each time step to ensure numerical stability. Finally, to account for the inherent variability in the chemical healing processes, a non-parametric probabilistic model is introduced and explicitly constructed onto the discrete operators (matrices) of the enthalpy evolution equation. Numerical examples illustrate the framework capability to simulate crack arrest and recovery. The stochastic analysis further quantifies the propagation of uncertainties from the chemical modeling level to the macroscopic healing efficiency, identifying the sensitivity of the structural recovery to variations in the healing reaction rates.