Predictive modelling and rheological investigation of titanium oxide–enhanced silica-based shear-thickening fluids

Abstract

A silica-based (SiO₂) shear-thickening fluid (STF) reinforced with titanium dioxide (TiO₂) nanoparticles was systematically examined across a broad range of temperature, shear rate, and applied stress. The TiO₂/SiO₂ system showed a distinct, stress-dependent, nonlinear shear-thickening behavior, with an optimal nanoparticle concentration that maximized viscosity increase. The synergistic interaction between SiO₂ and TiO₂ nanoparticles promoted the formation of dense hydroclaster structures, significantly enhancing resistance to flow at high shear rates and improving rheological stability across different thermal conditions. To accurately model and predict highly nonlinear viscosity behavior, advanced machine learning (ML) models, including Support Vector Regression (SVR), Random Forest (RF), XGBoost, and CatBoost, were developed and fine-tuned using Bayesian Optimization (BO). CatBoost outperformed others, achieving a coefficient of determination (R²) of 0.9991 and 0.9906 for training and testing datasets, respectively, along with a test root mean squared error (RMSE) of 7.2553 and a mean absolute error (MAE) of 4.0291. Additional validation confirmed the model’s generalization, with an R² of 0.9842. Moreover, Shapley additive explanation (SHAP) analysis was used to improve interpretability, revealing shear rate and shear stress as the most influential factors in viscosity changes. Meanwhile, temperature and TiO₂ concentration contributed positively within an optimal range. The close match between experimental results and ML predictions demonstrates the robustness of the proposed data-driven approach. This research highlights the potential of ML-based methods for the intelligent design and optimization of high-performance STFs for advanced impact mitigation, energy absorption, and vibration control applications.

Graphical Abstract

Introduction

Shear-thickening fluids (STFs) are a distinctive category of non-Newtonian materials that exhibit a sudden, significant increase in viscosity when subjected to applied shear or impact forces. This rapid transition from a fluid-like to a solid-like state is generally attributed to the formation of transition particle clusters, often referred to as hydroclasters, which result from intensified interparticle interactions under external loading conditions. Notably, STFs exhibit reversible behaviour, enabling them to revert to their original liquid state once the applied force is removed. This variable change has established STFs as essential materials for applications that require impact resistance and energy dissipation, such as earthquake-vibration mitigation systems, isolators, and dampers.^1–3

STFs generally consist of a dispersed phase alongside a dispersing medium.^4,5 A typical STF is composed of nano-sized silica (SiO₂) particles uniformly distributed within a liquid medium like polypropylene glycol or polyethylene glycol. These systems effectively absorb shock energy, making them appropriate for various types of dynamic situations and applications.^6–8 Recent advances have focused on enhancing the performance of STFs by introducing functional groups into the dispersed phase or by utilizing chemically modified particles. Several factors, including particle size, shape, mass fraction, dispersion medium, and temperature, significantly influence the rheological behaviour of STFs.^9–13

Research efforts have explored a wide array of dispersed-phase materials to optimize STF performance. For instance, Wagner and Egres demonstrated that suspensions that contain anisotropic particles like plates and rods exhibit enhanced shear thickening effects related to spherical particles due to their asymmetrical configuration. ¹⁴ Similarly, studies incorporating multiwalled carbon nanotubes and nanofibers of cellulose STFs have shown significant improvements in rheological properties, including greater shear-thickening and enhanced energy dissipation. ¹⁵ The relationship between dispersed particles and the dispersing medium has also been shown to play a pivotal role. Attalla and Osman underlined that polymer chains and functional units inside the medium greatly affect the STF’s performance, notably in terms of molecular interactions and structure stability. ¹⁶

Temperature sensitivity is another critical factor influencing STF behaviour. Wagner et al. experimentally demonstrated that the viscosity of STFs decreases with increasing temperature, driven by changes in particle interactions and fluid properties. ¹⁷ However, achieving optimal performance often requires balancing multiple constraints, as noted by Jiang et al. A new STF was developed by replacing SiO₂ with polystyrene-acrylic acid (PS-AA) nanoparticles (polystyrene-acrylic acid), leading to an additional substantial shear-thickening effect with increasing nanoparticle absorption. ¹⁸ Likewise, He et al. combined silica nanoparticles with cetyltrimethylammonium bromide (CTAB) to create an STF with a remarkable shear thickening effect due to enhanced nanoparticle surface interactions. ¹⁹ Wang et al. have demonstrated that STF characteristics can be further improved by combining hybrid materials, such as halloysite nanotubes, with silica nanoparticles. ²⁰ In a recent contribution, Khan et al. ²¹ investigated silica-based STFs reinforced with hexagonal boron nitride (h-BN) nanosheets, where enhanced rheological performance was primarily attributed to lubrication-dominated interactions, high thermal conductivity, and controlled hydrocluster formation. In contrast to such lubrication-based mechanisms, oxide-based nanoparticles such as titanium dioxide (TiO₂) offer a distinct advantage due to their higher surface energy, which promotes stronger interparticle friction and contact interactions. This influences shear thickening through friction-driven particle clustering, particle jamming, and temperature-dependent solvation layer effects, leading to fundamentally different rheological responses. Furthermore, TiO₂ exhibits distinct rheological behaviour influenced by particle size and shape^22–25 and provides better thermal stability and resistance to aggregation than silica, making it a promising option for hybrid STF formulations. Therefore, the present study aims to systematically investigate TiO₂-enhanced silica-based STFs and to establish a robust predictive framework to capture their nonlinear rheological behaviour.

The viscosity of STFs and nanofluids can be quantified utilizing rheometric instruments; however, such procedures are expensive, labour-intensive, and necessitate extensive repetitions to effectively characterize viscosity over wide parameter ranges. ²⁶ Conventional empirical models assume a continuous relationship between viscosity and shear rate, which limits their relevance for STFs that exhibit non-linear and discontinuous responses, including abrupt shear thickening. Only a limited number of studies have attempted to construct computational or analytical models capable of reproducing the complete viscosity profile, encompassing both the shear-thinning regimes and the discontinuous thickening zone. Galindo-Rosales et al. proposed an apparent viscosity expression that aligned with several experimental datasets, and a subsequent work extended this formulation to incorporate temperature and concentration effects. ²⁷ Shende et al. introduced an empirical relation based on free-volume theory that identifies threshold parameters and characterizes STF rheology using measurable quantities. ²⁸

Although these phenomenological models facilitate the generalization of experimental observations, they rely on simplifying assumptions and do not fully capture the physical mechanisms governing non-Newtonian behaviour. The need for numerous auxiliary parameters and repeated computations further restricts their applicability. Only a limited number of studies have explored predictive or hybrid modelling strategies that integrate rheological measurements with data-driven approaches to describe the nonlinear and discontinuous characteristics of STFs.^29,30 Furthermore despite the increasing incorporation of functional nanoparticles, such as TiO₂, to tailor STF performance, the specific influence of these additives on shear-thickening intensity, thermal stability, and the underlying interparticle interactions remains insufficiently characterized.

These limitations collectively define a critical research gap: the absence of a robust, interpretable framework that can simultaneously capture the nonlinear rheological response of nanoparticle-enhanced STFs, particularly those incorporating TiO₂, while overcoming the structural constraints of conventional phenomenological models. Although machine learning (ML) methods offer a promising alternative by providing improved predictive accuracy relative to empirical formulations,^31–33 their application in this domain has been largely exploratory and lacks mechanistic interpretability. For instance, Potnuri et al. applied support vector regression (SVR) to estimate product distributions in microwave-assisted co-pyrolysis demonstrating the model’s capability to incorporate catalyst characteristics, process conditions, and feedstock properties. ³⁴ Similarly Sun et al. employed SVR to predict crude-oil viscosity based on compositional and physicochemical parameters illustrating the utility of data-driven models for reservoir-fluid characterization under varying thermodynamic conditions. ³⁵ However analogous methodologies have yet to be systematically applied to predict the viscosity of TiO₂- and SiO₂-based STFs, while maintaining transparency and physical coherence.

The present study addresses these gaps by systematically investigating the rheological behaviour of STFs formulated with SiO₂ and TiO₂ particles across varying mass fractions. The experimental work encompasses shear-thickening response, temperature dependence, and viscoelastic properties through steady-state and dynamic rheological tests. To complement this experimental foundation, advanced ML models, including SVR, random forest (RF), extreme gradient boosting regression (XGBoost), and categorical boosting regression (CatBoost), are optimized using Bayesian optimization with Gaussian process surrogates (BO_GP) to predict STFs’ viscosity as a function of key experimental variables. Shapley additive explanations (SHAP) is further employed to interpret model outputs and quantify the contribution of each input feature, enabling a transparent, physically coherent understanding of the ML-based viscosity predictions. Collectively, these findings provide critical insights into integrating STF systems into composite materials for advanced impact-resistance and vibration-control applications.

Experimental methods

The selected materials

The study utilized SiO₂, which was purchased from Shanghai, China’s Aladdin Biochemical Technology Co., Ltd. ²² The SiO₂ particles supplied by Usolf Chemical Technology Cooperation Limited (Shandong Province, China) had an average particle size of approximately 15 nm. Their relative densities ranged from 2.338 to 2.655, and the measured pH values ranged from 3.88 to 4.77. Polypropylene glycol (PPG200) was a clear, stable fluid at ambient temperature. The hydroxyl group value ranged from 511 to 624 mg KOH/g, and TiO₂, with a molecular weight of approximately 80, was purchased from China Shanghai. ²²

Shear thickening fluid (STF) preparation

STF has been prepared in multiple phases to guarantee optimal rheological performance and dispersion. Before testing, TiO₂ and SiO₂ were heated in a vacuum oven at 70°C for 15 hours to remove any remaining moisture.

Through vigorous mixing, TiO₂ and SiO₂ particles were slowly added to the beaker in specific amounts. Throughout the procedure, the mixture was oscillated using an ultrasonic oscillator. Additionally, in a water bath, the beaker was kept at a constant temperature. Once all the dispersed-phase materials had been added, the mixture was stirred and oscillated using an ultrasonic oscillator until a homogeneous dispersion of the particles was achieved within the PPG. The TiO₂/SiO₂-STF was then maintained in a vacuum oven at 110°C. Once all bubbles in the system were eliminated, a steady SiO₂/TiO₂-STF system was produced. The TiO₂ to SiO₂ different mix ratios in four kinds of TiO₂/SiO₂-STF are displayed in Table 1.

Table 1.

Proportion of TiO₂ to SiO₂ in STF system.

STF (ingredient)	STF (ingredient percentage)
Silicon dioxide	15 + 0
Silicon dioxide + titanium dioxide	15 + 1
	15 + 3
	15 + 5

Rheological testing

Using a rheometer (TA Instruments AR2000, United States), the rheological characteristics of the produced TiO₂/SiO₂-STFs with different TiO₂-to-SiO₂ ratios were assessed in both steady-state and dynamic conditions. The experiments were performed at 19°C, a controlled ambient temperature. In this research, the cone-plate gap was maintained at 52 μm, and all measurements were conducted using the cone-plate geometry (diameter = 25 mm, cone angle = 2°). To prepare for testing, STF samples were evenly distributed on the rheometer’s bottom plate. After that, the top cone plate was lowered until the desired 52 μm gap was reached. The rheometer’s built-in temperature regulator was used to maintain temperature control throughout the trials. To assess the outcome of temperature on the STF samples’ performance, their rheological behaviour was examined over a range of 15 to 50°C. The shear rate was varied from 0.1 to 1000 s⁻¹ to perform the steady shear rheological testing. As a result, shear-thickening behaviour, including viscosity variation and the critical shear rate under different stress conditions, could be evaluated.

To perform dynamic rheological testing, oscillatory shear stress was applied to the STF samples over a range of 0.1-1000 Pa. These experiments provide insight into the STFs’ viscoelastic properties, including their ability to effectively store and dissipate energy under both low and high-stress conditions. The integration of steady-shear and dynamic rheological analysis provided a comprehensive understanding of the TiO₂/SiO₂-STFs’ viscoelasticity, shear-thickening behaviour, and temperature sensitivity. These results are essential for improving STF formulations for sophisticated uses.

Machine learning methods

Researchers have employed various ML techniques to predict the rheological and viscosity-related properties of complex fluids. In this study, the viscosity behavior of TiO₂/SiO₂-based STFs is modeled using both individual and ensemble learning approaches. SVR is applied as a particular method, while ensemble-based algorithms, RF, XGBoost, and CatBoost, are implemented to improve predictive performance through bagging and boosting frameworks. Hyperparameter tuning for all models was performed using BO_GP to ensure efficient and robust convergence. These models and the BO_GP optimization were selected for their proven effectiveness in capturing nonlinear, multivariate interactions and for their robustness in materials science prediction tasks. ³⁶

Support vector regressor

The Support Vector Machine (SVM), introduced by Corinna Cortes and Vladimir Vapnik in 1995, ³⁷ is a supervised learning algorithm that constructs an optimal separating hyperplane by maximizing the margin between classes. Its strong generalization capability, particularly in high-dimensional spaces, and its ability to handle nonlinear relationships through kernel functions make it well-suited for complex pattern recognition problems. ³⁸

In the regression context, referred to as SVR, the objective is to determine a function that approximates the target values while controlling model complexity. This is achieved by minimizing a regularized risk function under an ε-insensitive loss, where deviations within a predefined threshold ε are not penalized. The SVR model predicts a continuous output variable using input data represented by K-dimensional feature vectors $x_i$ and corresponding targets $y_i$ .

The SVR prediction function is expressed as equation (1):

$f\left(x\right)={\displaystyle \sum }_{i=1}^N\left({\alpha }_i-{\alpha }_i^{*}\right)\hspace{0.17em}K\left(x_i,x\right)+b$

(1)

where ${\alpha }_i$ and ${\alpha }_i^{*}$ are Lagrange multipliers, $x_i$ denotes the support vectors (i.e., training samples with nonzero multipliers), $K\left(x_i,x\right)$ is the kernel function that maps the input space into a higher-dimensional feature space, and $b$ is the bias term. ³⁹

Random forest

The Random Forest (RF) algorithm, an ensemble learning method based on Classification and Regression Trees (CART), was formally introduced by Leo Breiman in 2001. ⁴⁰ It has since been widely applied across diverse domains due to its strong predictive performance and robustness. As a statistical learning technique, RF employs bootstrap resampling (bagging) to generate multiple subsets from the original dataset. ⁴¹ Each subset is used to train an individual decision tree, with an additional layer of randomness introduced through feature subsampling at each split.

The final prediction is obtained by aggregating the outputs of all individual trees. For regression tasks, this aggregation is typically performed by averaging, expressed as equation (2):

$f\left(x^{\prime }\right)=\frac{1}{B}{\displaystyle \sum }_{b=1}^B{f}_b\left(x^{\prime }\right)$

(2)

where $f\left(x^{\prime }\right)$ denotes the predicted output for an unseen sample $x^{\prime }$ , $B$ is the total number of trees (or bootstrap samples), and $f_b\left(x^{\prime }\right)$ represents the prediction of the $b$ -th decision tree.

This ensemble strategy enhances model stability and generalization by reducing variance without substantially increasing bias. Contrary to common misconceptions, bagging primarily reduces variance, rather than increasing it, thereby mitigating overfitting. The use of both bootstrap sampling and random feature selection ensures that individual trees are decorrelated, which improves the robustness of the aggregated model. Consequently, RF is particularly effective for handling high-dimensional data and complex nonlinear relationships while being relatively insensitive to noise and outliers. ⁴¹

Extreme gradient boosting regression

XGBoost is a highly optimized, scalable implementation of the gradient boosting framework, specifically designed to improve predictive performance and computational efficiency. Unlike traditional gradient boosting algorithms that rely solely on first-order derivatives, XGBoost leverages both first- and second-order derivatives (gradients and Hessians) of the loss function, enabling more precise and stable parameter estimation at each boosting iteration. ⁴²

To control overfitting and enhance generalization, XGBoost introduces an explicit regularization term into its objective function. This includes both L1 (Lasso) and L2 (Ridge) penalties applied to the leaf weights of the decision trees. The regularization framework penalizes model complexity by constraining the number of leaves and the magnitude of leaf weights, thereby discouraging overly complex tree structures. Additionally, a gain-based splitting criterion is employed to evaluate candidate splits, ensuring an optimal balance between model fit and complexity.

From a computational perspective, XGBoost is engineered for efficiency and scalability. It supports parallel and distributed computation, utilizes cache-aware data structures to improve memory access patterns, and incorporates sparsity-aware algorithms to handle missing or sparse data efficiently. These optimizations enable significantly faster training and superior performance compared to traditional gradient boosting implementations, particularly for large-scale and high-dimensional datasets. ⁴³

The objective function of XGBoost is defined as equation (3):

$\mathcal{L}\left(\varphi \right)={\displaystyle \sum }_{i=1}^{n}l\left(y_i,{\widehat{y}}_i\right)+{\displaystyle \sum }_{k=1}^K\mathrm{\Omega }\left(f_k\right)$

(3)

where $l\left(y_i,{\widehat{y}}_i\right)$ is a differentiable loss function measuring the discrepancy between the observed value $y_i$ and the predicted value ${\widehat{y}}_i$ , and $\mathrm{\Omega }\left(f_k\right)$ is the regularization term for the $k$ -th tree.

The regularization term is given by as presented in equation (4):

$\mathrm{\Omega }\left(f\right)=\gamma T+\frac{1}{2}\lambda {\displaystyle \sum }_{j=1}^T{w}_j^2$

(4)

where $T$ denotes the number of leaves in a tree, $w_j$ represents the weight assigned to the leaf $j$ , $\gamma$ is a regularization parameter controlling the minimum loss reduction required to perform a split, and $\lambda$ is the L2 regularization parameter that penalizes large leaf weights. The L1 regularization component can also be incorporated to further promote sparsity in the model. ⁴⁴

CatBoost

In 2017, Yandex introduced its ML algorithm, CatBoost, which is based on Gradient Boosting Decision Trees (GBDT). ⁴⁵ CatBoost models operate by constructing a series of decision trees, each progressively reducing the loss. Configuration settings initially determine the tree generation rate, while an overfitting detector helps maintain model stability and prevent erratic behavior. ³⁷

The model updates in CatBoost are driven by an iterative objective function, mathematically represented as equation (5):

$h_t=\arg \underset{h\in \mathcal{H}}{\min }\mathbb{E}\left[\mathrm{\Phi }\left(y,F^{\left(t-1\right)}\left(x\right)+h\left(x\right)\right)\right]$

(5)

where $F^{\hspace{0.17em}\left(t-1\right)}\left(x\right)$ represents the accumulated output from previous iterations, $\mathrm{\Phi }\left(y,\widehat{y}\right)$ is the loss function, and $\mathbb{E}[·]$ denotes the expected value. The weak learner $h\left(x\right)$ is defined as presented in equation (6):

$h\left(x\right)={\displaystyle \sum }_{j=1}^T{w}_j\mathbf{1}\left(x\in Q_j\right)$

(6)

here, $w_j$ is the predicted value for the region $Q_j$ , $Q_j$ represents tree leaf regions, $\mathbf{1}(·)$ is the indicator function, and $T$ is the total number of leaves.

Building a tree in CatBoost involves several steps: computing candidate splits, converting categorical features to numerical representations, selecting the tree structure, and computing the leaf values. To improve prediction accuracy, CatBoost predominantly uses greedy algorithms. Features are organized by splits and assigned to each leaf in a consistent structure. Before constructing a new tree, feature permutations are applied rather than simple random shuffling. Pre-modelling parameters set the tree depth and constrain the selection of structure. When constructing subsequent trees, CatBoost evaluates performance based on a metric that reflects the improvement in the objective function. ⁴⁶

Bayesian optimization

Hyperparameter optimization is a critical step in ML, traditionally approached through methods like grid search and random search. ⁴⁷ Grid search systematically evaluates every possible combination of hyperparameters within the parameter space to find the optimal set, but it becomes computationally prohibitive as the parameter space expands due to its exponential complexity. ⁴⁸ Random search addresses this limitation by sampling a subset of parameter combinations, offering improved efficiency and suitability for high-dimensional datasets. However, both methods lack adaptability and rely on exhaustive or randomized search strategies.

Bayesian Optimization (BO) offers a more intelligent alternative by leveraging prior knowledge to guide the search. By combining the previous probability of the objective function with observed data, BO iteratively updates a posterior probability distribution to identify hyperparameters that minimize the objective function. This adaptive approach ensures that each iteration refines the search, focusing computational resources on promising regions of the parameter space. BO not only accelerates the finding of optimal hyperparameters but also reduces manual effort and improves model performance, making it a powerful tool for automating hyperparameter tuning in ML. ⁴⁹

BO aims to identify the optimal hyperparameters that minimize validation errors. Given a defined hyperparameter space and an objective function targeting the reduction of validation error, the optimization process within the BO framework is mathematically represented as equation (7):

$x^{*}=\arg \underset{x\in \mathcal{X}}{\min }f\left(x\right)$

(7)

where $x^{*}$ is the optimal hyperparameter combination, $\mathcal{X}$ is the parameter space, and $f\left(x\right)$ is the objective function. BO is a global optimization method for black-box functions, where the analytical form and derivatives are unknown. It leverages Bayes’ theorem to construct a probabilistic surrogate model of the objective function. ⁵⁰ BO then adaptively selects the next evaluation point based on prior observations, enabling a more efficient search for the optimal solution.

Shapley additive explanations

Tree-based ensemble learning methods like RF inherently offer feature relevance by counting how often a feature appears during model tree construction ⁴⁴ While evaluating feature importance through various methods reveals valuable insights, it doesn’t fully capture the overall impact of features or the complex relationships (both positive and negative) between features and the output. Therefore, understanding the full significance of the model’s features and their connections to the target variable is essential. This can be achieved through SHAP analysis, as proposed by Lundberg and Lee. ⁵¹

SHAP analysis uses Shapley values, derived from game-theoretic principles, to quantify each feature’s contribution. Global feature importance is typically obtained by averaging the absolute Shapley values across all instances in the dataset. These values are then visualized in descending order, where each point represents the Shapley value for a specific feature and instance. The x-axis represents the Shapley values (impact on model output), while the y-axis lists the features. The color gradient reflects the feature value (e.g., low to high), not the level of importance. Additionally, SHAP dependence plots capture the effects of feature interactions on the target variable, providing more detailed insights than traditional partial dependence plots. ⁵²

Lundberg and Lee developed several SHAP analysis methods, including DeepSHAP, LinearSHAP, KernelSHAP, and TreeSHAP. This study uses TreeSHAP. The explanation model for a prediction is formulated as equation (8):

$g\left(z^{\prime }\right)={\varphi }_0+{\displaystyle \sum }_{i=1}^M{\varphi }_i{z}_i^{\prime }$

(8)

where $g$ is the explanation model, $z^{\prime }$ represents the binary feature vector (presence or absence of features), $M$ is the total number of features, and ${\varphi }_i$ denotes the feature attribution values, with ${\varphi }_0$ as the base value.

Feature attribution for each feature is calculated as equations (9) and (10):

${\varphi }_i={\displaystyle \sum }_{S\subseteq N\backslash \left\{i\right\}}\frac{\mid S\mid !\left(M-\mid S\mid -1\right)!}{M!}\left[f\left(S\cup \left\{i\right\}\right)-f\left(S\right)\right]$

(9)

$f\left(S\right)=\mathbb{E}\left[f\left(x\right)\mid x_S\right]$

(10)

where $S$ is a subset of features, $N$ is the complete set of features, and $\mathbb{E}\left[f\left(x\right)\mid x_S\right]$ represents the expected model output conditioned on the subset $S$ .

Statistical Metrics

Three statistical metrics were employed to evaluate the efficacy of the ML models used in this investigation: coefficient of determination (R²), Mean Absolute Error (MAE), and Root Mean Square Error (RMSE). The following equations define these metrics, equations (11)–(13):

$R^2=1-\frac{{\sum }_i{\left(v_i^E-v_i^P\right)}^2}{{\displaystyle \sum }_i{\left(v_i^E-\overline{v}\right)}^2}$

(11)

$RMSE=\sqrt{\frac{1}{N}{{\sum }_i\left(v_i^E-v_i^P\right)}^2}$

(12)

$MAE=\frac{1}{N}{\sum }_i\left|v_i^E-v_i^P\right|$

(13)

where $v_{\mathrm{i}}^{\mathrm{E}}$ denotes the experimental viscosity for the ith data sample, while ${\mathrm{v}}_{\mathrm{i}}^{\mathrm{P}}$ represents the predicted viscosity for the same sample. N indicates the total number of data samples, and $\overline{v}$ signifies the average viscosity value across the dataset. For each statistical measure, the R² value was used to assess the difference between the predicted and actual data. At the same time, the MAE was calculated as the average absolute error between the predicted and actual values. The RMSE reflected the average of the errors between the projected and experimental data. ML models demonstrate better performance when they exhibit higher R² values and lower MAE and RMSE. For optimal predictive performance, a model’s R² value should approach 1.0, and both MAE and RMSE should ideally be 0.

Result and discussion

Steady-shear rheological analysis of STFs

The TiO₂/SiO₂-STF viscosity graphs stay within the 0.1 to 1000 s shear rate range during the steady-state test. Each sample was prepared at four different mass ratios, with PPG200 as the dispersive medium. The mass portions of TiO₂ used were 0%, 1%, 3%, and 5%, with the 15% STF made from SiO₂ serving as the baseline, as illustrated in Table 2. The viscosity curves for TiO₂/SiO₂-STF with various TiO₂ to SiO₂ ratios are displayed in Figure 1. This figure indicates that the SiO₂/TiO₂-STF system yielded a greater apparent value and an equally strong shear thickening effect in comparison to the STF, where the dispersed phase comprised just SiO₂. Incorporating dispersed-phase nano-TiO₂ improved the system’s overall viscosity, especially whenever the shear rate passed the point of critical shear rate. The TiO₂ nanoparticles may have gathered around the functional groups, preventing the system’s particle clusters from freely moving. Furthermore, the decreased mean distance between SiO₂ particles promoted the development of particle clusters, which raised the system’s viscosity.

Table 2.

Important coordinate for different TiO₂/SiO₂-STF. ²²

% Of SiO2/% of TiO2	Initial (viscosity/shear rate)	Critical (viscosity/shear rate)	Peak (viscosity/shear rate)
15% SiO2+ 0% TiO2	2.89 Pas 0.20 s−1	3.7 Pas 9 s−1	242.3 Pas 25.1 s−1
15% SiO2+ 1% TiO2	5.1 Pas 0.20 s−1	7 Pas 4.0 s−1	300 Pas 16.75 s−1
15% SiO2+ 3% TiO2	3.4 Pas 0.20 s−1	3.6 Pas 6.8 s−1	178.1 Pas 24.40 s−1
15% SiO2+ 5% TiO2	3.5 Pas 0.20 s−1	3.8 Pas 12.2 s−1	167.6 Pas 41.9 s−1

Figure 1.

Viscosity versus shear rate for the SiO₂/TiO₂-STF sample at numerous SiO₂/TiO₂ proportions. (a) Logarithmic Graph (b) Semi-Logarithmic Graph. ²²

When the proportion of dispersed TiO₂ particles rose, the system’s viscosity also rose; the curve was not constant, as demonstrated in Figure 1(a) and (b). For the TiO₂ mass fraction ratio, 1% is chosen as the optimal cut-off threshold. The SiO₂/TiO₂-STF system’s initial viscosity reduced, the apparent highest value rose, and the critical shear rate improved as the TiO₂ mass ratio increased whenever the TiO₂ percentage ratio passed this edge, which implies that a minor amount of TiO₂ boosted the silica particles’ interface, resulting in a greater shear-thickening impact and increased viscosity. Increasing the TiO2 proportion to 3% and then to 5% decreased the top viscosity and increased the critical shear rate, as observed with 1% TiO₂. This indicates that increasing the TiO₂ concentration in the fluid, encouraging fluid flow, and generating a lubricating effect reduce internal friction. This is because TiO₂ particles act as a dispersion agent, stabilizing the fluid’s silica particles and reducing their direct contact, thereby lowering viscosity and increasing the critical shear rate. These results suggest that the addition of TiO₂ has a complex effect on the system’s shear-thickening behaviour. Even though 1%TiO₂ provides the best viscosity and shear thickening characteristics, dispersion stabilization and lubrication effects may cause peak viscosity to decrease, and the shear thickening effect may not constantly be enhanced by increasing the TiO₂ concentration beyond this point.

Temperature sensitivity of STF

Figure 2(a) and (b) illustrates the correlation between the apparent viscosity and strain magnitude of the 1% SiO₂/STF at different test conditions temperatures. The following temperatures were used for the environmental testing: 15°C, 20°C, 30°C, 40°C, and 50°C. The shear rate’s scanning range is 0.1-1000 s⁻¹.

Figure 2.

TiO₂/SiO₂/STF sample viscosity versus shear rate over temperature. (a) Logarithmic Graph (b) Semi-Logarithmic Graph. ²²

The viscosity is similarly affected by the temperature shift, as seen in Figure 2(a) and (b) for the 1% STF system. Beginning and peak viscosities, fall, and the critical shear rate all increase with temperature, causing the viscosity curve to move downward accordingly. Because temperature rises have a considerable impact on the ST effect, it is noteworthy that, at every measured temperature, the system demonstrated exceptional stability with a 1% STF. At 15°C, the 1% STF system showed a roughly 36.20% reduction in the critical shear rate and a 27% increase in maximum viscosity compared with the controlled sample. At 50°C, the critical shear rate declined by up to 25.42%, while the viscosity of the 1% STF system significantly increased, reaching up to 51% higher than the 0% STF sample, as illustrated in Figure 2(a) and (b). However, this doesn’t lead to improved thermal stability. This is because temperature rises have a significant impact on the ST effect. In summary, the 1% STF consistently outperforms the 0% STF in terms of peak viscosity, critical shear rate, and resistance to temperature fluctuations. Adding only 1% TiO₂ to the shear-thickening fluid markedly improves its thermal stability, making it well-suited for applications that demand consistent performance in dynamic and elevated-temperature environments. These outcomes highlight that the TiO₂-based shear-thickening fluid has the potential to enhance STF performance and is a promising shear-thickening material.

ML model development and assessment

To ensure the robustness and generalization of the developed ML models, the complete dataset generated from the rheological experiments was first partitioned, with 15% of the observations reserved as an independent validation set. The remaining 85% of the data was then randomly split into a 70% training subset and a 15% testing subset. Each subset contained the same set of input variables and a single output variable (viscosity, η). The input features consisted of nanoparticle concentration (wt.%), shear rate (s⁻¹), temperature (°C), and the TiO₂ mass ratio, all of which govern the shear-thickening behaviour of the SiO₂/TiO₂-based system STFs. The input features were normalized to the (0, 1) range using min-max scaling for visual comparison, as shown in Figure 3, which displays the distributions of the training, validation, and test datasets. The consistency in distribution shapes across all three subsets demonstrates that the data partitioning was performed effectively, minimizing sampling bias and ensuring that the models are trained and tested on statistically similar populations.

Figure 3.

Distribution of independent features and target viscosity.

Figure 3 also presents the distribution of the target variable, viscosity. The histogram with kernel density estimation reveals a positively skewed distribution. This is expected for STFs, as most data points represent lower-viscosity states, with fewer observations capturing the high-viscosity shear-thickening regime. The violin-box plot confirms that the median viscosity is consistent across the training, validation, and test datasets. The elongated upper tails and outliers in this plot represent the critical, high-viscosity responses that are essential for accurately modelling the extreme shear-thickening behaviour. These distributions confirm that the dataset spans both low- and high-shear viscosity regimes, enabling the ML models to learn from the full rheological spectrum of the TiO₂/SiO₂-STFs.

Hyperparameter optimization using Bayesian optimization Gaussian process

In this study, BO_GP was implemented to optimize the hyperparameters of the SVR, RF, XGBoost, and CatBoost models. The optimization aimed to minimize the five-fold cross-validated RMSE across the training set, ensuring robustness and mitigating overfitting.

The search space for each algorithm was defined based on its influential hyperparameters, and the optimization process iteratively updated the Gaussian surrogate model to predict promising regions in the parameter space. At each iteration, the acquisition function (expected improvement) guided the selection of candidate parameter sets to balance exploration and exploitation. The final optimized parameters obtained for each model are presented in Table 3. The BO_GP convergence curves depicted in Figure 4 illustrate the progressive reduction in RMSE with increasing iterations. These convergence profiles validate the efficiency of BO_GP in navigating the hyperparameter space, ensuring that each model attained its optimal configuration for viscosity prediction.

Table 3.

Optimized parameters found by BO_GP for SVR, RF, XGBoost, and CatBoost.

Parameters	SVR RBF	Random forest	XGBoost	CatBoost
Regularization	C = 9433.7614	max_features = 0.7811	reg_alpha = 0.6785	l2_leaf_reg = 1.0446
		min_samples_split = 2	reg_lambda = 0.6147	min_data_in_leaf = 14
		min_samples_leaf = 1	gamma = 2.1763
			min_child_weight = 5
Tree structure	-	max_depth = 10	max_depth = 4	depth = 4
		n_estimators = 134	n_estimators = 227	iterations = 497
Learning	gamma = 0.0981	max_samples = 0.9905	learning_rate = 0.1788	learning_rate = 0.0716
	epsilon = 0.3442		subsample = 0.8256	subsample = 0.6805
			colsample_bytree = 0.7647	colsample_bylevel = 0.5905

Figure 4.

Bayesian optimization, Gaussian process five-fold CV convergence during training.

For the SVR model, the BO_GP optimization identified a high penalty parameter (C = 9433.76), indicating that the model required strong regularization to minimize training error. The kernel-specific parameters gamma (0.0981) and epsilon (0.3442) were tuned to control smoothness and margin tolerance, thereby ensuring accurate fitting of nonlinear rheological relationships.

In the RF model, the optimizer selected a maximum depth of 10 and 134 estimators, suggesting a moderately deep ensemble with sufficient trees for stability. The optimized values of max_features = 0.7811, min_samples_split = 2, and max_samples = 0.9905 indicate a balance between variance control and computational efficiency, enabling effective averaging across diverse decision trees.

The XGBoost model was optimized with reg_alpha = 0.6785 and reg_lambda = 0.6147, demonstrating the need for both L1 and L2 regularization to prevent overfitting. The learning rate (0.1788) and maximum depth (4) provided a good trade-off between model complexity and convergence stability, while subsample (0.8256) and colsample_bytree (0.7647) enhanced generalization through random feature and sample selection.

The CatBoost model achieved the lowest cross-validated RMSE via BO_GP optimization, as shown in Figure 4. The model converged after approximately 60 iterations, reaching a minimum RMSE of 7.937 and outperforming other algorithms. The optimal parameters included depth = 4, iterations = 497, learning_rate = 0.0716, and subsample = 0.6805, confirming the model’s ability to capture nonlinear dependencies effectively while preventing overfitting. Additionally, l2_leaf_reg = 1.0446 and colsample_bylevel = 0.5905 provided appropriate regularization and a balanced feature sampling, reinforcing CatBoost’s robustness in handling complex rheological patterns.

Prediction results

The predictive capabilities of the four BO_GP optimized SVR, RF, XGBoost, and CatBoost models were comprehensively evaluated across training, testing, and validation datasets. The performance indicators, including R², RMSE, and MAE, are summarized in Tables 4 and 5. Among all models, the CatBoost BO_GP model demonstrated the best performance, achieving R² values of 0.9991, 0.9906, and 0.9842 on the training, testing, and validation datasets, respectively. The corresponding RMSE (2.46-9.88) and MAE (1.51-5.03) values were the lowest among all models, confirming its superior predictive accuracy and generalization capability. The strong performance of CatBoost BO_GP can be attributed to its efficient ordered boosting algorithm, which minimizes overfitting while capturing complex nonlinear interactions among rheological parameters, including shear rate, nanoparticle concentration, TiO₂ ratio, and temperature.

Table 4.

Models’ performance during training and testing.

Model	Dataset	R2	RMSE	MAE
SVR BO_GP	Training	0.97	14.5290	5.5195
	Test	0.9515	16.4657	5.9714
RF BO_GP	Training	0.9975	4.2139	1.4381
	Test	0.9855	9.0143	3.3904
XGBoost BO_GP	Training	0.9997	1.4995	0.9269
	Test	0.9829	9.7755	4.5962
CatBoost BO_GP	Training	0.9991	2.4637	1.5106
	Test	0.9906	7.2553	4.0291

Table 5.

Models’ performance on validation data.

Model	R2	RMSE	MAE
SVR BO_GP	0.9692	13.8058	6.2643
RF BO_GP	0.9457	16.7601	5.9743
XGBoost BO_GP	0.9568	16.3507	7.5182
CatBoost BO_GP	0.9842	9.8787	5.0252

The XGBoost and RF models also exhibited high predictive accuracy, reflecting their effectiveness in modeling nonlinear systems. However, both showed slightly higher error values and marginally lower R² scores compared to CatBoost, suggesting reduced sensitivity to extreme viscosity variations. The SVR model, while performing reasonably well, showed comparatively higher errors, indicating its limited ability to capture abrupt rheological transitions characteristic of shear-thickening behavior.

The comparative scatter plots of predicted versus experimental viscosity values for all models across the training, testing, and validation datasets are presented in Figure 5. These plots provide a direct visualization of each model’s predictive fidelity and generalization capability. Ideally, a perfect prediction aligns along the 45° parity line, representing complete agreement between measured and predicted values. Among the models, the CatBoost BO_GP model exhibits the most coherent and compact distribution of data points along the parity line across all datasets, demonstrating exceptional consistency between predicted and experimental viscosities. The minimal deviation and tight clustering confirm the model’s strong generalization and its ability to accurately capture both low- and high-viscosity regimes.

Figure 5.

Scatter comparison of train, test, and validation for CatBoost, XGBoost, RF, and SVR.

The residual distributions illustrated in Figure 6 further validate these observations. The CatBoost BO_GP model displays a narrow, symmetric residual distribution centered around zero, indicating unbiased predictions and minimal systematic error. Most residuals fall within a range of −30 to 45, confirming the model’s stability and reliable performance across all datasets. The XGBoost BO_GP and RF BO_GP models exhibit wider residual spreads, with a slight tendency toward random deviations at higher viscosities. In contrast, the SVR BO_GP model exhibits the broadest and most uneven distribution, reflecting less stable predictive behaviour.

Figure 6.

Comparison of actual viscosity value and model predictions of train, test and validation, and residual histograms.

The Taylor diagram in Figure 7 provides a unified assessment of model performance by comparing correlation, standard deviation, and centered RMSD across the entire dataset. The CatBoost BO_GP model is positioned nearest to the reference point (R = 1, σ/σ_ref = 1), indicating the highest correlation and the most accurate replication of experimental viscosity variations. This close alignment confirms CatBoost’s superior capability to capture both the trend and amplitude of the rheological data with minimal deviation.

Figure 7.

Taylor diagrams of models on the entire dataset, CatBoost interpretations.

Prediction model interpretation

The SHAP analysis was performed to quantitatively interpret feature importance and directional effects in the BO_GP CatBoost model. The analysis was computed using the trained model on the testing set to evaluate the contribution of each input variable to the predicted viscosity of STFs.

As shown in Figure 8, the mean absolute SHAP values indicate that Stress (58.51) is the most influential feature, consistent with the literature on h-BN, which found that stress has a substantial influence on predicting STF viscosity ²¹ Stress substantially outweighs Shear rate (30.60), while Temperature (10.25) and TiO₂ concentration (4.01) demonstrate comparatively minor impacts. The beeswarm summary plot further elucidates the functional relationships: elevated Stress values correspond to strongly positive SHAP values, indicating a pronounced increase in predicted viscosity. In contrast, lower Stress values consistently reduce the model output. The relationship between Shear rate and viscosity is nonlinear, with both intermediate and high values exerting a negative influence on viscosity predictions. In contrast, Temperature and TiO₂ concentration exhibit limited effect magnitudes and broader distributions of values, indicating secondary roles, though their optimal ranges contribute positively to the output.

Figure 8.

Shap interpretation: Mean SHAP value and features distributions.

The SHAP dependence plots in Figure 9 elucidate the individual and interactive effects of features on the model’s viscosity predictions. Figure 9(a) reveals a strong positive relationship between Stress and the SHAP value. This relationship is moderated by Temperature, with the positive effect of Stress being most pronounced at lower temperature levels. Figure 9(b) shows a nonlinear relationship between shear rate and stress. SHAP values are positive at low to moderate shear rates but become negative beyond a critical threshold. This behavioural shift demonstrates a significant interaction with Stress, with combined high-stress and high-shear-rate conditions producing the most important adverse effect on viscosity predictions. Figure 9(c) indicates that increasing Temperature generally reduces the predicted viscosity. This negative correlation is particularly evident under low-stress conditions, indicating that Temperature acts as a moderating variable, attenuating stress-induced thickening. Figure 9(d) confirms that TiO₂ concentration has a non-monotonic relationship with viscosity. The maximum positive contribution occurs at approximately 1% concentration, consistent with experimental observations of optimal nanoparticle loading. Beyond this concentration, the reinforcing effect diminishes, especially under combined conditions of high stress and elevated temperature. Collectively, these dependence analyses verify that the model has captured the fundamental shear-thickening physics, including the stress-activated transition, shear-rate-dependent response, thermal suppression, and the optimal nanoparticle loading effect characteristic of STF systems.

Figure 9.

SHAP features dependence plots of the test set of CatBoost.

Conclusions and recommendations

This study presents a systematic investigation of the complex viscosity behavior of a TiO₂- SiO₂ hybrid STF. The rheological performance is examined as a function of nanoparticle composition and temperature, with particular emphasis on the influence of TiO₂ incorporation. The addition of a nominal fraction of TiO₂ is found to significantly enhance the shear-thickening response. Among the formulations evaluated, the hybrid system containing 1 wt% TiO₂ demonstrates an optimal balance of peak viscosity, critical shear rate, and thermal stability. Although temperature variations induce measurable changes in rheological behavior, this formulation retains a consistent shear-thickening characteristic, indicative of reliable thermal durability. The inclusion of 1 wt% TiO₂ in the STF matrix underscores its potential for applications requiring energy absorption, vibration damping, and protective performance under dynamic loading conditions.

Complementing the STF’s experimental characterization work, a rigorous ML framework was employed to model the complex viscosity behavior. The CatBoost model, optimized via BO, proved to be the most robust predictor, achieving superior accuracy and generalization across all datasets. Subsequent SHAP analysis provided critical model interpretability, confirming that the model’s predictions are physically consistent with rheological principles. Stress was identified as the dominant positive feature, while Shear rate exhibited a significant, nonlinear influence. The minor roles of Temperature and TiO₂ concentration align with the experimental observation of strong thermal adaptability and a defined optimal nanoparticle ratio.

The findings of this investigation confirm that the synergy between experimental characterization and ML modelling enables both the identification of optimal STF formulations and the development of accurate predictive tools. This integrated approach establishes a validated framework for advancing the design of functional fluids.

The analysis is based on a single dispersing medium, polyethylene glycol, and therefore does not account for the potential influence of alternative carrier fluids, such as silicon oil or ethylene glycol, on shear-thickening behavior. In addition, the predictive framework developed in this work is limited to rheological parameters derived from existing datasets, and application-oriented performance indicators such as impact absorption or energy dissipation efficiency were not included in the current modeling scope. Although the CatBoost model demonstrated superior performance among the tested algorithms, its predictions exhibited residual errors ranging from −30 to 45 Pa·s, indicating potential for further refinement. Future research should address these limitations by expanding the available dataset to encompass a broader range of nanoparticle concentrations, wider temperature intervals, and alternative dispersing media, thereby enabling the development of more comprehensive and transferable predictive models. Exploring advanced neural network architectures and physics-informed learning approaches could significantly improve predictive accuracy and extend model applicability to untested concentrations and conditions, ultimately enabling more efficient design of advanced STF systems.