Bulk-Surface Partial Differential Equations (BSPDEs) have become a hot topic over the last decade, as they arise in many applications across all fields of science. A BSPDE is a system where one or more PDEs posed in a compact domain $\Omega \subset R^d$ is/are coupled to one or more surface PDEs posed on the boundary $\partial \Omega$, which is assumed to be a smooth manifold. The theoretical analysis and numerical treatment of BSPDEs entail several open questions and research directions. On the other hand, the Virtual Element Method (VEM) is a recent polygonal counterpart of the Finite Element Method that has quickly gained popularity since the seminal work in [Beirao da Veiga et al., 2013, M3AS] thanks to its geometric flexibility. The first application of the VEM to BSPDEs was presented in [Frittelli et al., 2021, Numer Math], and the results were confined to the lowest order case. Here we present a high-order Virtual Element Method for BSPDEs in $d=3$ space dimensions. The proposed method relies on general polytopal meshes with both flat and curved faces. The flat faces are polygonal and with at most one curved edge. The curved faces, used to approximate the surface, are triangular. The method then couples an isoparametric FEM space on the surface with a novel bulk VEM space for the considered type of meshes. We illustrate our findings with numerical examples, including applications to battery modelling.