This is the most important step and great care and a good understanding of the technique are necessary if good results are required. This is particularly important for deep imaging as an error at the 0.01% level will significantly effect the photometry of the faintest sources.
There are two cases. In relatively blank fields, the sky is created from the object frames themselves. For crowded fields or large objects, the sky is created from frames that were specifically taken to measure the sky.
For deep exposures, the sky is computed from a subset of exposures and there will be one sky frame for each object frame. For accurate photometry, it is very important that the object frame is not included in the frames that are used to compute the sky. This is a weakness of the current eclipse recipe jitter. For H and K band observations the sky frame should be computed from frames that were taken 5-10 arc minutes either side of the object frame. For J band observations, these numbers can be doubled. For conditions where the sky background is varying rapidly (clouds or observations taken just after evening twilight) a more rapid sampling of the sky is necessary.
All sky frames contain objects, so one has to combine them with suitable clipping. A robust method is to first scale frames to a common median and then remove objects by rejecting the highest and lowest pixels. Rejecting the two highest and two lowest pixels would produce even more robust results. The remaining pixels are then averaged (the median can also be used, but it is a noisier statistic). The resulting sky frame is then scaled (to match the median of the object frame) and subtracted.
A more sophisticated approach is to do the sky subtraction in two steps. The first step reduces the data as described above, produces the combined image and then locates all the objects. In the second step, the data is rereduced with the knowledge of where the objects are. These objects are then excluded when the sky is estimated in the second pass. This is the approach used by the XDIMSUM package in IRAF and for very deep imaging it is the recommended package.