I am trying to run a linear regression model (ideally) to see whether age (continuous variable) affects levels of stress hormone (also continuous, dependent variable), i.e. hypothesis testing. My hormone level results are right-skewed and non-normal, and the residuals of a simple linear model have a similar distribution (see plots of residuals below). Age is normally distributed.
Scatter of raw data (has a couple outliers): Histogram of residuals:
I have also tried various transformations of my dependent variable (hormone levels) to improve the outcome of the model. The best normality plots came from the Tukey transformation. See plots below.
Scatter of transformed hormone level on age: Histogram of residuals:
However, the transformation yields a lower multiple R-squared and adjusted R-sqaured value than the original.
Further to this post (What if a transformed variable yields more normal and less heteroskedastic residuals but lower $R^2$?), I tried a model where the hormone level is log transformed and fitted with the exponent of age to correct for the log transformation. This is because I didn't know what the opposite of Tukey transformation would be and log seems to be the next best transformation for my data. The normality of residuals is good but the R2 values are even lower.
Given this information, which model is best - the one with best normality of residuals and lowest R2 or worst normality of residuals and highest R2?
Alternatively, is there another form of regression better suited to this data, given that it is continuous and has non-normal residuals?