0

I have a data frame 'dat' like below,

  structure(list(Tair = c(14.34, 14.14, 14.1, 13.92, 13.97, 14.1, 
  14.18, 14.15, 14.21, 14.38, 14.48, 14.49, 14.44, 14.11, 14.76, 
  12.46), Tsoil = c(15.59, 7.31, 7.31, 7.48, 7.83, 7.83, 7.83, 
  7.85, 7.87, 7.93, 7.96, 7.85, 7.79, 7.44, 7.47, 8.13), Treat2 = c("A", 
  "D 40%", "D 40%", "D 40%", "A", "A", "A", "D 40%", "D 40%", "D 40%", 
  "D 40%", "A", "D 50%", "A", "A", "D 50%"), LM.flux = c(1.29106069426062, 
   0.586392893610504, 0.541301472576533, 0.913764582996219, 1.28632065762852, 
   0.803288804475987, 1.1630166713827, 0.448368288220239, 0.438561131908755, 
   0.542394298341447, 0.520008006530466, 0.994437948111783, 0.727418651220693, 
   1.20284923105522, 0.867929974187552, 0.482064493861347), SWC = c(0.2397, 
   0.2408, 0.2408, 0.2705, 0.2408, 0.253, 0.2526, 0.2035, 0.2039, 
   0.233, 0.2326, 0.2382, 0.123, 0.1942, 0.2182, 0.2063)), row.names = c("51", 
   "214", "383", "552", "719", "887", "1057", "1227", "1396", "1568", 
   "1733", "1905", "2077", "2248", "2419", "2585"), class = "data.frame")


mode2_level <- paste("modrun2_SWC vs LM.flux", sep="")

modrun2 <- lm(LM.flux~SWC, data=dat, na.action=na.exclude, weights=(1/(LM.flux+0.49)),method = "qr")
summary(modrun2)

Now my point is how to find the relationship between SWC and LM.flux under different Treat2 treatments. Hope someone could help. Thanks

LEE
  • 316
  • 2
  • 8
  • 1
    Perhaps this would help - [Linear Regression and group by in R](https://stackoverflow.com/questions/1169539/linear-regression-and-group-by-in-r) – arg0naut91 Sep 26 '20 at 21:02
  • 1
    Hi @LEE, i think there are many good answers in the link provided by arg0naut91.Please give them a try. The answer provided below, it should work too. – StupidWolf Sep 26 '20 at 21:24
  • 1
    Thanks all of you, that link helped a lot. – LEE Sep 27 '20 at 07:14

1 Answers1

0

Something like this?

library(dplyr)
library(tidyr)

mod_df <- 
  dat %>% 
  group_by(Treat2) %>% 
  nest() %>% 
  rowwise() %>% 
  mutate(mod = list(lm(LM.flux ~ SWC, data = data, na.action = na.exclude,
                       weights=(1/(LM.flux+0.49)), method = "qr"))) %>% 
  ungroup()

mod_df
#> # A tibble: 3 x 3
#>   Treat2 data             mod   
#>   <chr>  <list>           <list>
#> 1 A      <tibble [7 x 4]> <lm>  
#> 2 D 40%  <tibble [7 x 4]> <lm>  
#> 3 D 50%  <tibble [2 x 4]> <lm>

mod_df$mod
#> [[1]]
#> 
#> Call:
#> lm(formula = LM.flux ~ SWC, data = data, weights = (1/(LM.flux + 
#>     0.49)), na.action = na.exclude, method = "qr")
#> 
#> Coefficients:
#> (Intercept)          SWC  
#>       1.391       -1.393  
#> 
#> 
#> [[2]]
#> 
#> Call:
#> lm(formula = LM.flux ~ SWC, data = data, weights = (1/(LM.flux + 
#>     0.49)), na.action = na.exclude, method = "qr")
#> 
#> Coefficients:
#> (Intercept)          SWC  
#>     -0.7744       5.7764  
#> 
#> 
#> [[3]]
#> 
#> Call:
#> lm(formula = LM.flux ~ SWC, data = data, weights = (1/(LM.flux + 
#>     0.49)), na.action = na.exclude, method = "qr")
#> 
#> Coefficients:
#> (Intercept)          SWC  
#>       1.090       -2.945
Iaroslav Domin
  • 2,698
  • 10
  • 19
  • @laroslav Domin why all the coefficients are the same after trying this code (mod_df$mod)? – LEE Sep 28 '20 at 09:01
  • @LEE what do you mean by saying "all the coefficients are the same". `mod_df$mod` shows three different linear models with different coefficients as can be seen in my answer. – Iaroslav Domin Sep 29 '20 at 18:04