An h-adaptive scheme of first-order least-squares meshfree method is presented. A posteriori error estimates, which can be readily computed from the residual, are also presented. For elliptic problem the error indicators are further improved by applying the Aubin-Nitsche method. In the proposed refinement scheme, Voronoi cells are utilized to insert nodes at appropriate positions. Through numerical examples, it is demonstrated that the error indicators reveal good correlations with the actual errors and the adaptive first-order least-squares meshfree method is effectively applied to the localized problems such as the shock formation in fluid dynamics.