This paper is concerned with the construction of optimized sparse grid approximation spaces for elliptic pseudodifferential operators of arbitrary order. Based on the framework of tensor-product biorthogonal wavelet bases and stable subspace splittings, we construct operator-adapted subspaces with a dimension smaller than that of the standard full grid spaces but which have the same approximation order as the standard full grid spaces, provided that certain additional regularity assumptions on the solution are fulfilled. Specifically for operators of positive order, their dimension is