Probabilistic approach to analysis of synergism in mathematical models of biochemical networks is introduced. It is based on system analysis concept when information on the importance of a parameter of a complex biochemical model is evaluated as part of joint interaction with a complete set of model parameters. For example, this approach accounts for uncertainties in the estimates of enzyme activities and kinetic parameters involved in kinetic modelling of the networks and/or concentration of metabolite or cofactors involved in the interaction of a pathway with perturbations on a cellular level. The parameters are considered as random variables with assumed corresponding probability distribution functions, and total effects of their variability on the network fluxes are evaluated. A numerical measure of synergism of an individual parameter with respect to interaction with model parameters is defined as the difference between the ensemble expected value of conditional variance for the complementary parameters and the variance of conditional expected value of the particular parameter relative to the total parameter ensemble dispersion. In order to demonstrate the concept, the proposed method is applied to two simple cases and to a complex model. The first case is the analysis of synergism between activator and substrate in uni-uni type I mechanism. In the second example, synergism between enzymes involved in flux through a serial pathway is evaluated. As an example of a complex system, synergism between glycogenolytic flux in a skeletal muscle and involved cellular level cofactors is analyzed.