correlation with ordinal variable

Question

I am trying to analyze experiments where small groups of animals are all treated with a drug at different doses. So each value on the x-axis represents one group and each value on the y axis the measurements in a single animal.

I want to figure out if the drug treatment has an effect on the values I am measuring. I already ran ANOVA and Dunnet Tests to generate the 'between group' stats. However, now I want to use all groups together to see if there is an overall relationship between dose and response. So I believe this means that treatment is now an ordered categorical variable (ordinal).

The two questions I have are:

1. How should I decide if I should treat the drug dose as a common numeric variable (based on the exact dose given) rather than a categorical (or ordered numeric) variable?
2. How do I interpret the different results between running a linear model with the ordinal variable (none, low, med, high) vs turning it into a ordered numeric (1,2,3,4).

Below is a quick rundown of the data and how I am analyzing it so far. I have a related question on different correlation stats link that I will post independently so this post stays a little focused. The data is provided at end of post.

#// overview
library(tidyverse)
ggplot(data = stats) + geom_boxplot(mapping = aes(x = ord, y = value)) +
                       geom_point(mapping = aes(x = ord, y = value))

To see if there is a correlation between drug dose and the measured value, I can turn the categorical variable into ordered integers and then simply run a linear model:

summary(lm(value ~ card, stats))
#> 
#> Call:
#> lm(formula = value ~ card, data = stats)
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -459.36 -194.56  -14.33  131.75  603.84 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)   1714.7      108.5  15.801  < 2e-16 ***
#> card          -291.4       38.9  -7.491 8.98e-09 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 272.8 on 35 degrees of freedom
#> Multiple R-squared:  0.6159, Adjusted R-squared:  0.6049 
#> F-statistic: 56.11 on 1 and 35 DF,  p-value: 8.983e-09

ggplot(data = stats) + geom_smooth(mapping = aes(x = card, y = value), method = "lm") +
                       geom_point(mapping = aes(x = card, y = value))
#> `geom_smooth()` using formula 'y ~ x'

Or I can use the drug dose (low = 150, med = 500, high = 1500) and assume that there is a linear correlation between dose given and exposure in the animal and then run a linear fit.

summary(lm(value ~ num, stats))

#> 
#> Call:
#> lm(formula = value ~ num, data = stats)
#> 
#> Residuals:
#>    Min     1Q Median     3Q    Max 
#> -353.3 -214.1  -20.6  153.4  729.5 
#> 
#> Coefficients:
#>               Estimate Std. Error t value Pr(>|t|)    
#> (Intercept) 1297.63661   61.02984  21.262  < 2e-16 ***
#> num           -0.56801    0.07401  -7.674 5.26e-09 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 268.7 on 35 degrees of freedom
#> Multiple R-squared:  0.6273, Adjusted R-squared:  0.6166 
#> F-statistic:  58.9 on 1 and 35 DF,  p-value: 5.261e-09

ggplot(data = stats) + geom_smooth(mapping = aes(x = num, y = value), method = "lm") +
                       geom_point(mapping = aes(x = num, y = value))
#> `geom_smooth()` using formula 'y ~ x'

However, it turns out that I can also directly run a linear model on the ordinal data, but the output is very different and I don't really know how to interpret the results:

summary(lm(value ~ ord, stats))
#> 
#> Call:
#> lm(formula = value ~ ord, data = stats)
#> 
#> Residuals:
#>     Min      1Q  Median      3Q     Max 
#> -416.34 -200.00  -39.24  139.01  646.86 
#> 
#> Coefficients:
#>             Estimate Std. Error t value Pr(>|t|)    
#> (Intercept)  1380.32      85.02  16.234  < 2e-16 ***
#> ord2_low     -239.50     132.50  -1.808  0.07981 .  
#> ord3_medium  -423.05     120.24  -3.518  0.00129 ** 
#> ord4_high    -910.91     120.24  -7.576 1.03e-08 ***
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 268.9 on 33 degrees of freedom
#> Multiple R-squared:  0.6482, Adjusted R-squared:  0.6162 
#> F-statistic: 20.26 on 3 and 33 DF,  p-value: 1.254e-07

#// data using dput()

stats <- structure(list(ord = c("1_none", "1_none", "1_none", "1_none", 
"1_none", "1_none", "1_none", "1_none", "1_none", "1_none", "2_low", 
"2_low", "2_low", "2_low", "2_low", "2_low", "2_low", "3_medium", 
"3_medium", "3_medium", "3_medium", "3_medium", "3_medium", "3_medium", 
"3_medium", "3_medium", "3_medium", "4_high", "4_high", "4_high", 
"4_high", "4_high", "4_high", "4_high", "4_high", "4_high", "4_high"
), num = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 150, 150, 150, 150, 
150, 150, 150, 500, 500, 500, 500, 500, 500, 500, 500, 500, 500, 
1500, 1500, 1500, 1500, 1500, 1500, 1500, 1500, 1500, 1500), 
card = c(1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2, 
2, 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 
4, 4), value = c(1472.43, 1083.59, 1277.04, 2027.18, 1761.27, 
1027.95, 1309.47, 1879.03, 963.98, 1001.27, 876.03, 859.17, 
975.37, 937.42, 1233.27, 1569.03, 1535.48, 1001.83, 871.19, 
826.28, 1221.87, 1165.72, 972.36, 944.66, 715.28, 1096.28, 
757.28, 411.13, 267.38, 622.18, 599.03, 531.28, 568.38, 287.04, 
665.28, 312.27, 430.17)), row.names = c(NA, -37L), class = c("tbl_df", 
"tbl", "data.frame"))

^{Created on 2021-01-17 by the reprex package (v0.3.0)}

Answer 1

With respect to the last issue, you treat them as comparisons. R has made the "none" the base category, so we interpret it as follows:

For animals with no treatment, as compared to those with low treatment, there is a -239.50 change in value, on average. Repeat the comparison with None vs other condition. You can change what becomes the base "comparison" category by setting the variable as a factor and specifying the base.