Unsaturated soil has some resistance for slope failure, but this resistance decreases due to water infiltration and seepage. The material point method (MPM) has made remarkable progress for slope failure analysis. However, MPMs for unsaturated remain relatively few, and most existing approaches have employed explicit time integration, which suffers from high computational costs and pore water pressure instability. To solve these problems, we have recently developed a semi-implicit MPM with fractional-step method for unsaturated soil [1]. Nevertheless, there is still some room for improvement. One of them is suppressing a numerical error caused by using linear shape functions. While B-spline shape functions mitigate this error, in the vicinity of the physical boundary of simulation target, only a part of the shape function support domain is often covered empty domain and is not determined uniquely. As a result, the nodes near the boundaries of the spatially fixed grid do not coincide with the physical boundaries, leading to numerical errors. In order to resolve such inconveniences, this study aim to enhance the MPM for unsaturated soil by modified B-spline shape functions near the boundary using weighted least squares (WLS) [2]. Nodes are classified into inside, boundary and outside nodes based on the number of particles per grid, and the shape functions of boundary nodes are modified by replacing those of empty nodes with a WLS corrected function. We implement them in the semi-implicit MPM to stabilize the calculation of deformation analysis of unsaturated soil. To demonstrate the capability and performance of the proposed method, and validation analysis is carried out for simulating a model experiment of seepage failure.