twoway lfit out of y-axis range

[Daisy Lee]

Dear All,

Hello. I am creating a scatter plot and lfit by using the below code. Then I get the plot below the code.

twoway (scatter y x_normalized, msymbol(oh) (lfit y x_normalized), legend(off) scheme(white_tableau) ylabel(0(0.2)1)

As you can see, the lfit line goes below the y-axis value of 0, so the y-axis doesn’t meet the x-axis. I want the lfit line to be plotted only between the y-axis values of 0 and 1.

This dataset is US county-level. Variable y is not affected by the county population size. But the variable x is seriously affected by the population size and has a large variation. (min = 0, max = 825.4105)
Therefore, I used variable 'x_normalized' which is defined by x/population in the above code and plot. I think I fully addressed the problem of population size, but the y-axis problem in the plot occurs whether I use x or x_normalized.

How can I address the problem of lfit going out of y-axis range?

Thank you for your help.

p.s. If I omit the observations with the top 10% x value or assign weights based on the county population, the problem does not occur.

[Nick Cox]

The issue here is statistical rather than graphical. A linear fit is not a good idea for these data. You could truncate the line by using twoway function to show the line for a specified range of x but a much better solution would be to fit a more suitable curve such as a generalized linear model with logarithmic link.

[Nick Cox]

My previous was written when my only internet access was using my phone.

Here's a positive suggestion. You hide what y and x are except something for US counties. But it seems clear that y is by definition positive or just possibly positive or zero. In that context an exponential fit that never predicts a positive mean is immensely more plausible than a linear fit that will go negative at some point.

glm y x_normalized, link(log) f(poisson) vce(robust) 

predict fitted 

scatter y x_normalized || mspline fitted x_normalized, bands(200) legend(off)

is a sketch of code to try.

If y is also bounded by 1, then

​​​​​​​glm y x_normalized, link(logit) f(binomial) vce(robust)

is more nearly correct in principle, but is likely to give similar results for your data, as most outcome values are much nearer 0 than 1, and in that circumstance logit and log are close.