The paper gives a detailed error analysis of the boundary element method applied to the solution of the first biharmonic boundary value problem (1) \(\Delta^ 2\phi =0\) in \(\Omega \subset {\mathbb{R}}^ 2\); \(\phi =f_ 1\) and \(\partial_ n\phi =f_ 2\) on \(\Gamma\). The boundary equations are derived on the basis of the direct integral equation method proposed earlier by S. Christiansen and P. Hougaard [J. Inst. Math. Appl. 22, 15-27 (1978; Zbl 0391.35028)]. The error analysis is made for spline approximations defined by the Galerkin approach as well as by the collocation technique. Although the system of boundary equations is not strongly elliptic, optimal order error estimates for the unknown boundary densities are obtained in some Sobolev spaces via a perturbation technique. Using some integral representation formula, the authors derive pointwise error estimates of the form \(| \phi (x)-\phi_ h(x)| \leq ch^ p\) for the solution of (1).