We have exploited a variety of techniques to study the universality and stability of the scaling properties of Harper's equation, the equation for a particle moving on a tight-binding square lattice in the presence of a gauge field, when coupling to next-nearest sites is added. We find, from numerical and analytical studies, that the scaling behavior of the total width of the spectrum and the multifractal nature of the spectrum are unchanged, provided the next-nearest-neighbor coupling terms are below a certain threshold value. The full square symmetry of the Hamiltonian is not required for criticality, but the square diagonals should remain as reflection lines. A bicritical line is found at the boundary between the region in which the nearest-neighbor terms dominate and the region in which the next-nearest-neighbor terms dominate. On the bicritical line a different critical exponent for the width of the spectrum and different multifractal behavior are found. In the region in which the next-nearest-neighbor terms dominate, the behavior is still critical if the Hamiltonian is invariant under reflection in the directions parallel to the sides of the square, but a new length scale enters, and the behavior is no longer universal but shows strongly oscillatory behavior. For a flux per unit cell equal to 1/q the measure of the spectrum is proportional to 1/q in this region, but if it is a ratio of Fibonacci numbers the measure decreases with a rather higher inverse power of the denominator.