对415个文件进行线性回归，只输出文件名、回归系数、显着性

Question

我是 R 的初学者，我正在学习分析一些生物数据的基础知识。我有 415 个 .csv 文件，每个文件都是一种真菌。每个文件有 5 列 - (YEAR, FFD, LFD, MEAN, RANGE)

    YEAR    FFD LFD RAN MEAN
1   1950    NA  NA  NA  NA
2   1951    NA  NA  NA  NA
3   1952    NA  NA  NA  NA
4   1953    NA  NA  NA  NA
5   1954    NA  NA  NA  NA
6   1955    NA  NA  NA  NA
7   1956    NA  NA  NA  NA
8   1957    NA  NA  NA  NA
9   1958    NA  NA  NA  NA
10  1959    140 141 1   140
11  1960    NA  NA  NA  NA
12  1961    NA  NA  NA  NA
13  1962    NA  NA  NA  NA
14  1963    NA  NA  NA  NA
15  1964    NA  NA  NA  NA
16  1965    155 156 1   155
17  1966    NA  NA  NA  NA
18  1967    NA  NA  NA  NA
19  1968    152 153 1   152
20  1969    NA  NA  NA  NA
21  1970    NA  NA  NA  NA
22  1971    161 162 1   161
23  1972    NA  NA  NA  NA
24  1973    143 144 1   143
25  1974    NA  NA  NA  NA
26  1975    NA  NA  NA  NA
27  1976    NA  NA  NA  NA
28  1977    NA  NA  NA  NA
29  1978    NA  NA  NA  NA
30  1979    NA  NA  NA  NA
31  1980    NA  NA  NA  NA
32  1981    NA  NA  NA  NA
33  1982    155 156 1   155
34  1983    NA  NA  NA  NA
35  1984    NA  NA  NA  NA
36  1985    157 158 1   157
37  1986    170 310 140 240
38  1987    173 274 101 232
39  1988    192 236 44  214
40  1989    234 320 86  277
41  1990    172 287 115 213
42  1991    148 287 139 205
43  1992    140 278 138 206
44  1993    152 273 121 216
45  1994    142 319 177 228
46  1995    261 318 57  287
47  1996    247 315 68  285
48  1997    164 270 106 230
49  1998    186 187 1   186
50  1999    235 236 1   235
51  2000    NA  NA  NA  NA
52  2001    309 310 1   309
53  2002    203 308 105 256
54  2003    140 238 98  189
55  2004    204 313 109 267
56  2005    253 313 60  287
57  2006    247 300 53  279
58  2007    185 295 110 225
59  2008    259 260 1   259
60  2009    296 315 19  309
61  2010    230 303 73  275
62  2011    247 248 1   247
63  2012    206 207 1   206
64  2013    NA  NA  NA  NA
65  2014    250 317 67  271

首先我想查看每个文件的回归系数（线的斜率），以及所有文件的显着性（p 值）。

我可以单独做：

fruit<-read.csv(file.choose(),header=TRUE)
yr<-fruit[,1]
ffd<-fruit[,2]
res<-lm(ffd~yr)
summary(res)

当我对数据执行此操作时，我得到：

Call:
lm(formula = ffd ~ yr)

Residuals:
    Min      1Q  Median      3Q     Max 
-77.358 -20.858  -5.714  22.494  96.015 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -4162.0710   950.1439  -4.380 0.000119 ***
yr              2.1864     0.4765   4.588 6.55e-05 ***
---
Signif. codes:  
0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 38.75 on 32 degrees of freedom
  (31 observations deleted due to missingness)
Multiple R-squared:  0.3968,    Adjusted R-squared:  0.378 
F-statistic: 21.05 on 1 and 32 DF,  p-value: 6.549e-05

目前我需要的唯一信息是回归系数 (2.1864) 和 p 值 (6.549e-05)

如果我可以让 R 循环浏览 415 个文件，并以包含 3 列的表格形式给出输出：文件名、回归系数和显着性，那么完美的输出将是。将有 415 行，每个文件一个。

然后我想做 YEAR~LFD、YEAR~RANGE 和 YEAR~MEAN。我希望我可以轻松地编辑 YEAR~FFD 的代码并为其他 3 个回归运行它。

Answer 1

以下代码可能会起作用。
我已经在两个文件中使用您的数据对其进行了测试。完成所有工作的函数如下：

regrFun <- function(DF){
  fit <- lm(DF[[1]] ~ DF[[2]])
  coef(summary(fit))[2, c(1, 4)]
}
regrList <- function(iv, L){
  res <- lapply(seq_along(L), function(i){
    dftmp <- L[[i]]
    cfs <- regrFun(dftmp[c(1, iv)])
    data.frame(file = names(L)[i], Estimate = cfs[1], p.value = cfs[2])
  })
  res <- do.call(rbind, res)
  row.names(res) <- NULL
  res
}

现在读入数据文件。在下面的代码行中，在明显的位置用通用文件名部分替换 "pattern"。

filenames <- list.files(pattern = "pattern")
df_list <- lapply(filenames, read.csv)
names(df_list) <- filenames

然后计算你想要的值。

results_list <- lapply(2:ncol(df_list[[1]]), regrList, df_list)
names(results_list) <- names(df_list[[1]][-1])

Answer 2

首先我模拟 5 个 csv 文件，其中的列看起来像你的：

for(i in 1:5){
  tab=data.frame(
      YEAR=1950:2014,
      FFD= rpois(65,100),
      LFD= rnorm(65,100,10),
      RAN= rnbinom(65,mu=100,size=1),
      MEAN = runif(65,min=50,max=150)
  )
  write.csv(tab,paste0("data",i,".csv"))
}

现在，我们需要一个包含您目录中所有文件的向量，这对您来说会有所不同，但请尝试使用模式参数以某种方式创建它：

csvfiles = dir(pattern="data[0-9]*.csv$")

所以我们使用了 tidyverse 的三个库，我猜每个 csv 文件并没有那么大，所以下面的代码读取所有文件，按源分组并执行回归，注意你可以使用调用列数据框，不必重命名：

library(dplyr)
library(purrr)
library(broom)

csvfiles %>% 
map_df(function(i){df = read.csv(i);df$data = i;df}) %>%  
group_by(data) %>% 
do(tidy(lm(FFD ~ YEAR,data=.))) %>% 
filter(term!="(Intercept)")

# A tibble: 5 x 6
# Groups:   data [5]
  data      term  estimate std.error statistic p.value
  <chr>     <chr>    <dbl>     <dbl>     <dbl>   <dbl>
1 data1.csv YEAR   -0.0228    0.0731    -0.311 0.756  
2 data2.csv YEAR   -0.139     0.0573    -2.42  0.0182 
3 data3.csv YEAR   -0.175     0.0650    -2.70  0.00901
4 data4.csv YEAR   -0.0478    0.0628    -0.762 0.449  
5 data5.csv YEAR    0.0204    0.0648     0.315 0.754  

您只需更改 lm(FFD ~ YEAR,data=.) 中的公式即可获得其他回归

Answer 3

data.table 版本，使用 StupidWolf 的 csv 文件布局和名称，具有请求的字段：

library(data.table)

input.dir = "/home/user/Desktop/My Folder/" # adjust to your needs
csvfiles <- list.files(path=input.dir, full.names=TRUE, pattern=".*data(.*)\\.csv") # adjust pattern

上面，我使用了更具体的正则表达式模式，但如果您想处理该文件夹中的所有 csv 文件，您可以选择 pattern="*.csv"。

# order the files
csvfiles <- csvfiles[order(as.numeric(gsub(".*data(.*)\\.csv", "\\1", csvfiles)))]

# function to read file and return requested columns
regrFun <- function(x){
    DT <- fread(x)
    fit <- lm(FFD ~ YEAR, data=DT)
    return(as.list(c(filename=basename(x), coef(summary(fit))[2, c(1, 4)])))
}

# apply function and rename columns
DT <- rbindlist(lapply(csvfiles, regrFun))
setnames(DT, c("filename", "regression coefficient", "significance"))

DT

结果：

        filename regression coefficient       significance
  1:   data1.csv     -0.113286713286712 0.0874762832713643
  2:   data2.csv     -0.044449300699302  0.457096760642717
  3:   data3.csv     0.0464597902097902  0.499618510612891
  4:   data4.csv     -0.032473776223776  0.638494798460044
  5:   data5.csv     0.0562062937062939  0.452955919860998
 ---                                                      
411: data411.csv     0.0381555944055959  0.544185411150829
412: data412.csv    -0.0672202797202807  0.314346452751388
413: data413.csv      0.116564685314687 0.0694785724198052
414: data414.csv    -0.0908216783216786  0.110811677724832
415: data415.csv    -0.0282779720279721  0.638766712090455

Answer 4

您可以编写一个在单个文件上运行的 R 脚本，然后通过终端在每个文件上运行它。

脚本只是一个包含代码的 .R 文件。 要在您将在终端上执行的每个文件上运行它（使用 bash）

for file in $(ls yourDataDirectory); do
    Rscript yourScriptFile.R $file >> finalOutput
done

这将在 yourDataDirectory 中的每个文件上运行 yourScriptFile.R 中的脚本，并将输出保存在 finalOutput 中。

脚本代码本身将与您已经编写的非常相似，但您将使用命令行传递的参数代替file.choose()，如here 所述，并且您只需要打印信息你很感兴趣，而不是summary 的输出。 如果您正确格式化脚本输出，finalOutput 甚至可以是 csv 文件。