Fundamental studies have been carried out experimentally and theoretically on the magnesiothermic reduction of silica with different Mg/SiO2molar ratios (1–4) in the temperature range of 1073 to 1373 K with different reaction times (10–240 min). Due to the kinetic barriers occurring in metallothermic reductions, the equilibrium relations calculated by the well-known thermochemical software FactSage (version 8.2) and its databanks are not adequate to describe the experimental observations. The unreacted silica core encapsulated by the reduction products can be found in some parts of laboratory samples. However, other parts of samples show that the metallothermic reduction disappears almost completely. Some quartz particles are broken into fine pieces and form many tiny cracks. Magnesium reactants are able to infiltrate the core of silica particles via tiny fracture pathways, thereby enabling the reaction to occur almost completely. The traditional unreacted core model is thus inadequate to represent such complicated reaction schemes. In the present work, an attempt is made to apply a machine learning approach using hybrid datasets in order to describe complex magnesiothermic reductions. In addition to the experimental laboratory data, equilibrium relations calculated by the thermochemical database are also introduced as boundary conditions for the magnesiothermic reductions, assuming a sufficiently long reaction time. The physics-informed Gaussian process machine (GPM) is then developed and used to describe hybrid data, given its advantages when describing small datasets. A composite kernel for the GPM is specifically developed to mitigate the overfitting problems commonly encountered when using generic kernels. Training the physics-informed Gaussian process machine (GPM) with the hybrid dataset results in a regression score of 0.9665. The trained GPM is thus used to predict the effects of Mg-SiO2mixtures, temperatures, and reaction times on the products of a magnesiothermic reduction, that have not been covered by experiments. Additional experimental validation indicates that the GPM works well for the interpolates of the observations.