Stochastic differential equations are widely used in various fields; in particular, the usefulness of duality relations has been demonstrated in some models such as population models and Brownian momentum processes. In this study, a discussion based on combinatorics is made and applied to calculate the expectation values of functions in systems in which evolution is governed by stochastic differential equations. Starting with the duality theory of stochastic processes, some modifications to the interpretation and usage of time-ordering operators naturally lead to discussions on combinatorics. For demonstration, the first and second moments of the Ornstein–Uhlenbeck process are re-derived from the discussion on combinatorics. Furthermore, two numerical methods for practical applications are proposed. One method is based on a conventional exponential expansion and the Padé approximation. The other uses a resolvent of a time-evolution operator, along with the application of the Aitken series acceleration method. Both methods yield reasonable approximations. Particularly, the resolvent and Aitken acceleration show satisfactory results. These findings will provide a new way of calculating expectations numerically and directly without using time-discretization.