# correlation with ordinal variable

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)