Differentiation matrices are an important tool in the implementation of the spectral collocation method to solve various types of problems involving differential operators. Fractional differentiation of Jacobi orthogonal polynomials can be expressed explicitly through Jacobi–Jacobi transformations between two indexes. In the current paper, an algorithm is presented to construct a fractional differentiation matrix with a matrix representation for Riemann–Liouville, Caputo and Riesz derivatives, which makes the computation stable and efficient. Applications of the fractional differentiation matrix with the spectral collocation method to various problems, including fractional eigenvalue problems and fractional ordinary and partial differential equations, are presented to show the effectiveness of the presented method.