We present a statistical mechanical theory for polymer-solvent systems based on integral equations derived from the polymer Kirkwood hierarchy. Integral equations for pair monomer-monomer, monomer-solvent, and solvent-solvent correlation functions yield polymer-solvent distribution, chain conformation in three dimensions, and scaling properties associated with polymer swell and collapse in athermal, good, and poor solvents. Variation of polymer properties with solvent density and solvent quality is evaluated for chains having up to 100 bonds. In good solvents, the scaling exponent upsilon has a constant value of about 0.61 at different solvent densities computed. For the athermal solvent case, the gyration radius and scaling exponent decrease with solvent density. In a poor solvent, the chain size scales as N-upsilon with the value of the exponent being about 0.3, compared with the mean field value of 1/3.