Linearly Regressing Pooled Data

I'm trying to regress some simple pooled data. My data has 60 observations and three columns: Weight, Height, and Sex (female=1, male=0).

If I regress thus, Weight ~ Height + Sex, my model is fairly satisfactory, but the residuals are not homoscedastic (green errors are male, blue female):

I tried regressing on the log of Weight and/or Height, but that didn't do much. What should I do to make the residuals homescedastic and/or make my model more accurate? Any help would be appreciated.

Edit

Doing a generalized regression model gives the following.

Generalized least squares fit by REML
Model: Weight ~ h + s
Data: P149
AIC      BIC    logLik
514.2221 524.4374 -252.1111

Variance function:
Structure: Different standard deviations per stratum
Formula: ~1 | Sex
Parameter estimates:
0         1
1.0000000 0.6685307

Coefficients:
Value Std.Error   t-value p-value
(Intercept)  27.197499  51.88129  0.524226  0.6022
h             1.852382   0.75634  2.449128  0.0174
s           -25.284478   5.53300 -4.569755  0.0000

Correlation:
(Intr) h
h -0.997
s -0.524  0.466

Standardized residuals:
Min         Q1        Med         Q3        Max
-1.6655243 -0.6879858 -0.1839396  0.5628971  3.9857544

Residual standard error: 22.13369
Degrees of freedom: 60 total; 57 residual


With this s. residual plot:

Could someone please explain how precisely this model is different from a standard multiple regression model? Thanks.

• You should really consider adding an interaction term to the model. I would expect the relationship between height and weight to be different between the sexes. So, you should at least test that. – Roland Feb 19 '14 at 7:59

You should use a model framework that allows you to model the different variation in males and females. Package nlme provides a generalized least squares implementation (gls) that can be used in conjunction with the varIdent function to specify the appropriate variance structure. E.g., you could try something like
library(nlme)

• Thanks. This is the R output. Is gls correcting the variance or just telling me the structure? Generalized least squares... Model: Weight ~ Height + Sex AIC BIC logLik 514.2221 524.4374 -252.1111 Variance function: Structure: Different standard deviations per stratum Formula: ~1 | Sex Parameter estimates: 0 1 1.0000000 0.6685307 Coefficients: Value Std.Error t-value p-value (Intercept) 27.197499 51.88129 0.524226 0.6022 Height 1.852382 0.75634 2.449128 0.0174 Sex -25.284478 5.53300 -4.569755 0.0000 – N4v Feb 18 '14 at 14:26