We study the Fredholm integro-differential equation $$D^{2s}_xσ(x) + \int_{-∞}^{+∞} k(x y)σ(y)dy = g(x)$$ by the wavelet method. Here $σ(x)$ is the unknown fun tion to be found, $k(y)$ is a convolution kernel and $g(x)$ is a given function. Following the idea in [7], the equation is discretized with respect to two different wavelet bases. We then have two different linear systems. One of them is a Toeplitz-Hankel system of the form $(H_n + T_n)x = b$ where $T_n$ is a Toeplitz matrix and $H_n$ is a Hankel matrix. The other one is a system $(B_n + C_n)y = d$ with condition number $k = O(1)$ after a diagonal scaling. By using the preconditioned conjugate gradient (PCG) method with the fast wavelet transform (FWT) and the fast iterative Toeplitz solver, we can solve the systems in $O(n$ ${\rm log}$ $n)$ operations.