/
experimental-design-advanced.qmd
1049 lines (805 loc) · 67.2 KB
/
experimental-design-advanced.qmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
```{r echo = FALSE}
pacman::p_load(tidyverse, readxl, knitr, kableExtra, Hmisc,
grid, agricolae, patchwork, desplot, modelr, see)
```
# Fortgeschrittene Designs {#sec-experimental-design-advanced}
*Letzte Änderung am `r format(fs::file_info("experimental-design-advanced.qmd")$modification_time, '%d. %B %Y um %H:%M:%S')`*
> *"The best time to plan an experiment is after you've done it." --- Roland Fisher*
In diesem Kapitel wollen wir uns mit der Erstellung von komplexeren experimentellen Designs beschäftigen. Wie immer nicht mit der Auswertung, aber ich habe dir Modellierungsvorschläge immer wieder angegeben.
- [Randomized complete block design (RCBD, 3-faktoriell)](#sec-rcbd-3fac)
- [Split plot design (3-faktoriell) oder Spaltanlage](#sec-split)
- [Subsampling](#sec-subsampling)
- [Incomplete block design (3-faktoriell und 2-faktoriell)](#sec-incomplete)
- [Strip plot design (3-faktoriell) oder Streifenanlage](#sec-strip-plot)
Dabei ist zu Bedenken, dass diese Designs sehr viel Aufwand in der Anlegen der Anlage sowie der praktischen Auswertung auf dem Feld bedürfen. Hier musst du dich dann nochmal beraten lassen, bevor du mit einem Experiment startest. Wir schauen uns jetzt aber erstmal einen Überblick über die gängigsten komplexeren experimentellen Designs an. Beachte auch die kurze Veröffentlichung über [Versuchsergebnisse interpretieren](https://www.dlg.org/de/landwirtschaft/themen/pflanzenbau/pflanzenernaehrung/dlg-kompakt-3-2021) der DLG e.V. - Deutsche Landwirtschafts-Gesellschaft. Dort findest du nochmal mehr Kontext über komplexere experimentelle Designs und Grundsätze für die Interpretation von Versuchsergebnissen, die du natürlich auch auf deine Planungen übertragen kannst. Bevor du die Ergebnisse interpretieren kannst musst du ja eine solide Versuchsplanung durchgeführt haben.
::: callout-important
## Auf ein Wort zum Umfang des Kapitels
Ich werde mich in diesem Kapitel für das Erste auf die wichtigen experimentellen Design beschränken. Wenn sich aus der Lehre oder Beratungstätigkeit noch Bedarf an weiteren experimentellen Designs ergibt, werde ich die Designs entsprechend ergänzen. Insbesondere hier gilt, dass [DSFAIR - Data Science for Agriculture in R](https://schmidtpaul.github.io/DSFAIR/DesigningExperiments.html) einen sehr umfangreichen Überblick über experimentelle Designs und deren Auswertung liefert.
**Wenn dir ein experimentelles Design fehlt, dann schreibe mir gerne eine Mail und ich schaue, dass ich das experimentelle Design ergänze.**
Mehr Informationen findest du auch durch das Inhaltsverzeichnis von dem R Paket `agricolae` mit [Experimental Designs with {agricolae}](https://myaseen208.com/agricolae/articles/ExperimentalDesign.html). Auch ist es sehr hilfreich sich Datenbeispiele auf [{agridat} - Datensätze mit Abbildungen in `{desplot}`](https://kwstat.github.io/agridat/reference/index.html) anzuschauen. Am meisten hilft vermutlich auch die Referenzliste des [R Paketes `{FielDHub}`](https://didiermurillof.github.io/FielDHub/reference/index.html). Leider sind die generischen Abbildungen in `{FielDHub}` nicht so toll, wie die Abbildungen, die du dir selber in `{desplot}` baust.
:::
Das charmante an dem neuen R Paket `{FielDHub}` ist die Verbindung einer Shiny App zusammen mit dem normalen Code, den du hier siehst. Ich habe dir immer die Referenz zu `{FielDHub}` hinzugefügt, in der ist auch die Anleitung zu finden für die Shiny App. Um `{FielDHub}` im Browser zu starten, musst du nur folgenden Code in R ausführen, nachdem du das R Paket installiert hast.
```{r}
#| eval: false
library(FielDHub)
run_app()
```
In der @fig-exp-adv-1 sehen wir einmal die Übersicht über die drei häufigsten komplexeren experimentellen Designs mit dem RCBD (3-faktoriell), dem Split plot Design (3-faktoriell) sowie dem Subsampling eines RCBD (2-faktoriell). Die drei experimentellen Designs werden immer mal wieder etwas durcheinander gewirbelt. Du kannst dir im Prinzip die Sachlage wie folgt vorstellen. Im RCBD (3-faktoriell) werden die beiden Behandlungsfaktoren über die Blöcke randomsiert. Die beiden Behandlungen sind vollkommen durcheinander in der @fig-exp-adv-1-1. Im Fall eines Split plot Design (3-faktoriell) wird eine Behandlung innerhalb der anderen Behandlung über die Blöcke randomisiert. Wir haben also eine Spaltanlage vorliegen, da wir immer eine Behandlung in einer Spalte randomisieren. Das ist schon in der @fig-exp-adv-1-2 zu sehen. Die Farben bilden die zweite Behandlung ab. Innerhalb einer Spalte ist dann die erste Behandlung randomisiert. Das Subsampling nimmt jetzt eine Behandlung raus. Wir haben dann nur noch eine Behandlung als Faktor, die wir in den Blöcken randomisieren. Die Idee des Subsamplings ist nun, dass wir wiederholt Pflanzen innerhalb der Behandlung messen. Wir mitteln jetzt aber nicht über die Pflanzen in einem Block sondern nehmen alle Pflanzen mit in die statistische Auswertung. Das ist eigentlich die zentrale Idee des Subsamplings, wie wir das machen, ist dann die andere Frage. In der @fig-exp-adv-1-3 sehen wir schön wie sich die Pflanzen-ID's über das Feld verteilen. Die Farben stellen die vier Behandlungslevel dar.
```{r}
#| echo: false
#| warning: false
#| label: fig-exp-adv-1
#| fig-align: center
#| fig-height: 6
#| fig-width: 6
#| fig-cap: "Vergleich der drei häufigsten komplexeren experimentellen Designs mit dem RCBD (3-faktoriell), dem Split plot Design (3-faktoriell) sowie dem Subsampling eines RCBD (2-faktoriell). Die Zahlen bei dem Subsampling stellt die ID's der individuellen Beobachtungen dar."
#| fig-subcap:
#| - "Sechs Eisendüngerbehandlung (Text: fe_1 bis fe_6) sowie vier Bodenbehandlung (Farbe: schwarz, rot, orange, grün) in vier Blöcken."
#| - "Sechs Eisendüngerbehandlung (Text: fe_1 bis fe_6) sowie vier Bodenbehandlung (Farbe: schwarz, rot, orange, grün) in vier Blöcken."
#| - "Vier Bodenbehandlung (Farbe: schwarz, rot, orange, grün) mit den individuellen Beobachtungen in vier Blöcken."
#| layout-nrow: 1
#| column: page
trt_fac1_soil <- str_c("soil_", 1:4)
n_trt_fac1_soil <- n_distinct(trt_fac1_soil)
trt_fac2_fert <- str_c("fe_", 1:6)
n_trt_fac2_fert <- n_distinct(trt_fac2_fert)
n_block <- 4
## rcbd
fac2rcbd_out <- design.ab(trt = c(n_trt_fac2_fert, n_trt_fac1_soil),
design = "rcbd",
r = n_block,
seed = 42)
fac2rcbd_out$bookRowCol <- fac2rcbd_out$book |>
bind_cols(expand.grid(rows = 1:n_trt_fac2_fert,
cols = 1:(n_trt_fac1_soil*n_block))) |>
mutate(trt_fac2_fert = paste0("fe_", A),
trt_fac1_soil = paste0("soil_", B),
block = paste0("Block ", block))
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = trt_fac2_fert, cex = 1, shorten = "no", col=trt_fac1_soil,
data = fac2rcbd_out$bookRowCol,
main = "Randomized complete block design (3-faktoriell)",
show.key = FALSE, key.cex = 0.5)
## split plot
splitplot_out <- design.split(trt1 = trt_fac1_soil,
trt2 = trt_fac2_fert,
r = n_block,
seed = 42)
splitplot_out$bookRowCol <- splitplot_out$book |>
mutate(block = paste0("Block ", block),
cols = plots |> str_sub(2,3) |> as.integer(),
rows = splots |> as.integer())
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = trt_fac2_fert, cex = 1, shorten = "no", col=trt_fac1_soil,
data = splitplot_out$bookRowCol ,
main = "Split plot design (3-faktoriell)",
show.key = FALSE, key.cex = 0.5)
## subsampling
rcbd_long_tbl <- expand_grid(block = 1:4,
trt = 1:4,
rep = 1:6) |>
mutate(trt = factor(trt, labels = trt_fac1_soil),
block = factor(block, labels = str_c("Block ", 1:4)),
pid = 1:n(),
pid_text = str_pad(1:96, width = 2, pad = "0")) |>
group_by(block) |>
mutate(trt = sample(trt))
rcbd_plot_tbl <- rcbd_long_tbl |>
bind_cols(expand_grid(cols = 1:16, rows = 1:6))
desplot(data = rcbd_plot_tbl, flip = TRUE,
form = block ~ cols + rows | block,
out1 = trt, out1.gpar = list(col = "grey", lty = 1),
out2 = pid, out2.gpar = list(col = "grey", lty = 1),
main = "Subsampling (2-faktoriell)",
text = pid_text,
cex = 1, show.key = FALSE, shorten = "no", col = trt)
```
## Genutzte R Pakete
Wir wollen folgende R Pakete in diesem Kapitel nutzen.
```{r echo = TRUE}
#| message: false
pacman::p_load(tidyverse, magrittr, agricolae, dae, desplot,
janitor, FielDHub,
conflicted)
conflicts_prefer(dplyr::filter)
conflicts_prefer(magrittr::set_names)
cbbPalette <- c("#000000", "#E69F00", "#56B4E9", "#009E73",
"#F0E442", "#0072B2", "#D55E00", "#CC79A7")
```
An der Seite des Kapitels findest du den Link *Quellcode anzeigen*, über den du Zugang zum gesamten R-Code dieses Kapitels erhältst.
## Randomized complete block design (RCBD, 3-faktoriell) {#sec-rcbd-3fac}
Das Maximum des machbaren ist eigentlich ein dreifaktorielles Modell im *randomized complete block design*. Natürlich geht auch noch mehr, aber der Grundsatz ist eigentlich, dass wir uns maximal zwei Behandlungsfaktoren und dann ein bis zwei Cluster anschauen. Die Cluster sind dann meist einmal der klassische Block plus ein Faktor für verschiedene Lokalisationen. Wir nehmen hier jetzt einmal zwei Behandlungsfaktoren und dann noch einen klassischen Block dazu. Beide Behandlungsfaktoren sind ineinander und dann natürlich im Block genestet. Hier aber erstmal das Modell mit den drei Faktoren für den besseren Überblick. Wir haben einmal die Düngung `fert`, dann den Faktor Boden `soil` sowie die verschiedenen Blöcke durch `block`.
$$
y \sim \overbrace{soil}^{f_1} + \underbrace{fert}_{f_2} + \overbrace{block}^{f_3}
$$
mit
- dem Faktor Dünger `fert` und den sechs Düngestufen als Level `fe_1` bis `fe_6`.
- dem Faktor Boden `soil` und den vier Bodenarten als Level `soil_1` bis `soil_4`.
- dem Faktor Block `block` und den zwei Leveln `block_1` und `block_2`.
Die Struktur der Daten ist wie folgt gegeben. Jedes Level der Düngerstufe `fert` ist in dem Faktor des Bodens `soil` enthalten. Die beiden Behandlungen sind dann wiederum jeweils vollständig in den Blöcken `block` vorhanden. In der @fig-mermaid-adv-rcbd sehen wir einmal den Zusammenhang zwischen den drei Faktoren
```{mermaid}
%%| label: fig-mermaid-adv-rcbd
%%| fig-width: 6
%%| fig-cap: "text"
flowchart LR
A(fert):::factor --- B(((nestet))) --> E(block):::factor
C(soil):::factor --- D(((nestet))) --> E
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
```
Wir wollen uns jetzt die Daten einmal selber erstellen und das R Paket `agricolae` nutzen um uns ein passendes Design zu bauen. Abschließend schauen wir uns den Versuchsplan einmal mit `desplot` an.
::: callout-tip
## Modell zur Auswertung
Wir rechnen ein multiples lineares Modell für die statistische Analyse.
```{r}
#| eval: FALSE
fit <- lm(drymatter ~ soil + fert + soil:fert + block, data = rcbd_3f_tbl)
```
:::
### ... mit `expand_grid()`
Bei dem *randomized complete block design* haben wir ja Glück, wir müssen nur darauf achten, das die einzelnen Beobachtungen innerhalb der Blöcke vollständig randomisiert sind. Daher Randomisieren wir ganz am Ende einmal die Saplte `fert` und die Spalte `soil` durch. Wir hätten uns das `mutate()` auch sparen können und stattdessen einfach die Funktion `slice_sample(prop = 1)` nutzen. Die Funktion permutiert dann alles innerhalb der gruppierten Blöcke durch. Mach wie es dir besser gefällt und du es besser nachvollziehen kannst.
```{r}
three_fct_long_tbl <- expand_grid(block = 1:2, soil = 1:4, fert = 1:6) |>
mutate(block = factor(block, labels = str_c("Block ", 1:2)),
soil = factor(soil, labels = str_c("soil_", 1:4)),
fert = factor(fert, label = str_c("fe_", 1:6))) |>
group_by(block) |>
mutate(fert = sample(fert), # Randomisierung fert
soil = sample(soil)) # Randomisierung soil
three_fct_long_tbl
```
Jetzt mache ich es etwas anders bei der Erstellung des Grids der Spalten und Zeilen. Ich will später beim `desplot()` für die Blöcke durch `| block` separieren. Daher setzte ich die Spaltennummerierung `cols` jeweils auf 1 bis 4 und bilde die Blöcke durch ein davor geschaltetes `rep` ab. Wenn du die Variable `rep` nicht möchtest, kannst du auch die `cols` auf 1 bis 8 setzten, darfst dann aber nicht für die Blöcke durch `rep` separieren.
```{r}
three_fct_plot_tbl <- three_fct_long_tbl |>
bind_cols(expand_grid(rep = 1:2, cols = 1:4, rows = 1:6))
```
In der @fig-exp-adv-6 sehen wir dann einmal das *randomized complete block design* mit drei Faktoren einmal dargestellt. Wir sehen gut, wie die zwei Behandlungen vollständig randomisiert wurden. Beachte, dass die Farben den Faktor `soil` darstellen und die Labels dann den Faktor `fert`. Du kannst mit der Option `show.key = FALSE` auch die Legende ausschalten. Bei sehr komplexen Designs mit vielen Faktorstufen ist es dann doch mal ratsam, sich die Legende mit ausgeben zu lassen.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-adv-6
#| fig-align: center
#| fig-height: 5
#| fig-width: 7
#| fig-cap: "Randomized complete block design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` vollständig randomisiert in den beiden Blöcken."
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = fert, cex = 1, shorten = "no",
col = soil,
data = three_fct_plot_tbl ,
main = "Randomized complete block design (3-faktoriell)",
show.key = TRUE)
```
Auch hier können wir dann den Versuchsplan relativ einfach raus schreiben. Ich entferne noch die Spalte `rep`, da ich die Spalte nicht weiter brauchen werde. Dann wären wir auch schon mit dem Design fertig.
```{r}
#| eval: false
#| message: false
#| warning: false
three_fct_plot_tbl |>
select(-rep) |>
write_xlsx("template_sheet.xlsx")
```
### ... mit `{FielDHub}`
Für die Erstellung des *randomized complete block design* in der Shiny App schaue dir die Anleitung unter [Generates a Full Factorial Design mit `{FielDHub}`](https://didiermurillof.github.io/FielDHub/reference/full_factorial.html) einmal an. Hier einmal der Code für die Erstellung in R. Wir haben jetzt vier Behandlungen für den Faktor `soil` sowie sechs Behandlungen für den Faktor `fe` vorliegen und wollen jeweils drei Blöcke anlegen. Jetzt einmal die generische Funktion `full_factorial()` in der wir aber keine Namen der Level übergeben können. Ohne die Namen wird die Abbildung maximal wirr. Ich erkenne dann gar nichts mehr, wenn ich nur noch Zahlen stehen habe.
```{r}
rcbd_fac3_obj <- full_factorial(setfactors = c(4, 6),
reps = 3, l = 1, type = 2)
```
Deshalb gibt es die Möglichkeit sich ein Grid zu erstellen und anhand der Namen des Grid dann die Parzellen zu benennen. Leider ist das Vorgehen etwas suboptimal was das Einfärben gleich angeht.
```{r}
data_factorial <- lst(fct = c(rep(c("soil", "fe"), c(4, 6))),
lvl = c(str_c("soil_", 1:4), str_c("fe_", 1:6))) |>
as.data.frame()
data_factorial
```
Wir bauen uns jetzt die Daten mit der Funktion `full_factorial()` setzen aber die Option `setfactors = NULL` damit die Funktion die Faktoren und deren Level aus unserem Datensatz `data_factorial` nimmt. Das ist jetzt auch nur so semi praktisch, aber immerhin könnten wir hier beliebig viele Faktoren prinzipiell nehmen.
```{r}
rcbd_fac3_obj <- full_factorial(setfactors = NULL, reps = 3, l = 1, type = 2,
data = data_factorial)
```
Dann können wir uns einmal das erstellte Feldbuch wiedergeben lassen, was du dann auch in Excel raus schreiben kannst. Nicht immer brauchst du alle Informationen. Zum Beispiel ist der Ort `location`, wo das Experiment durchgeführt wird eher zweitrangig. Du könntest dir auch zusätzliche Daten generieren lassen, aber das kannst du dir dann auch selber mit Werten ausdenken.
```{r}
rcbd_fac3_obj$fieldBook |>
as_tibble()
```
Dann erstellen wir auch schon die @fig-exp-basic-fac3-fieldhub mit dem *randomized complete block design*. Leider ist die Abbildung nur noch an der Grenze zu okay. Wir haben wilde Farben, die jeweils die Faktorkombination wiedergeben. Leider helfen die Farben wenig und auch ist alles sehr unübersichtlich. Das Paket `{desplot}` kann definitiv mehr. Hier würde ich dann wirklich auf den Plot hier aus `{FielDHub}` verzichten und stattdessen mir den Plot wie oben in der @fig-exp-adv-6 aus den den Feldbuch selber bauen.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-basic-fac3-fieldhub
#| fig-align: center
#| fig-height: 4
#| fig-width: 10
#| fig-cap: "Schema der Pflanzungen nach den Behandlungen für das *randomized complete block design*. Die Farbwhal und die Beschriftung ist mehr als gewöhnungsbedürftig. Die Abbildung ist nicht sehr übersichtlich und nur begrenzt als Darstellung geeignet."
rcbd_fac3_obj |> plot()
```
### ... mit `{agricolae}`
Jetzt bauen wir unser Design nochmal in `{agricolae}` nach. Hier brauchen wir dann aber alles als Vektor, sonst wird es zu unübersichtlich in den Funktionen. Darüber hinaus ist es dann auch mal was anders und du siehst nochmal eine andere Art den Code zu schreiben. Prinzipiell hätten wir auch im vorherigen Teil alles erstmal in Vektoren lagern können. Aber gut, hier erstmal alle Level der Faktoren vorbereiten und die Anzahl der Blöcke auf vier gesetzt. Dann brauchen wir noch die Länge der Level, also die Anzahl an Düngerstufen und Bodenarten. Dafür nutzen wir die Funktion `distinct_n()`.
```{r}
trt_fac1_soil <- str_c("soil_", 1:4)
n_trt_fac1_soil <- n_distinct(trt_fac1_soil)
trt_fac2_fert <- str_c("fe_", 1:6)
n_trt_fac2_fert <- n_distinct(trt_fac2_fert)
n_block <- 4
```
Dann können wir schon das *randomized complete block design* mit drei Faktoren über die Funktion `design.ab()` erstellen. Wie immer heißt die Wiederholung auch hier `r`, das ist zwar sehr einheitlich aber manchmal auch verwirrend. Auch musst die beiden Vektoren aneinanderkleben damit die Funktion funktioniert.
```{r}
rcbd_fac3_obj <- design.ab(trt = c(n_trt_fac2_fert, n_trt_fac1_soil),
design = "rcbd",
r = n_block,
seed = 42)
```
Dann können wir auch schon unser Ergebnis der Funktion einmal aufarbeiten. Wie immer fehlt das Positionsgrid, so dass wir hier nochmal tätig werden müssen. Dann wollen wir noch die Faktoren wieder umbenennen, so dass auch die Level passen. Am Ende wähle ich noch die Spalten aus, die wir dann später brauchen werden. Meistens müssen wir auch bei den komplexeren Designs und der Nutzung von `agricolae` im Nachgang sehr viel selber machen. Teilweise lohnt es sich da für mich nicht, für zwei Zeilen Code eine Funktion zu nutzen, der ich dann auch noch sehr viel Nachprogrammieren muss. Man merkt hier eben auch das Alter von dem R Paket `agricolae`.
```{r}
rcbd_fac3_book_tbl <- rcbd_fac3_obj$book |>
bind_cols(expand.grid(rows = 1:n_trt_fac2_fert,
cols = 1:(n_trt_fac1_soil*n_block))) |>
mutate(trt_fac2_fert = str_c("fe_", A),
trt_fac1_soil = str_c("soil_", B),
block = paste0("Block ", block)) |>
select(block, fert = trt_fac2_fert, soil = trt_fac1_soil, rows, cols)
```
Natürlich hat die Funktion auch keinen `sketch`, wenn man ihn braucht. Das haben wir ja jetzt schon selber mit dem Positionsgrid programmiert. In der @fig-exp-adv-2 siehst du dann das Ergebnis des Versuchsplans. Wir haben hier nochmal zwei Blöcke zusätzlich zu dem obigen Beispiel genommen, dann siehst du nochmal schöner, wie sich die beiden Behandlungen in den Blöcken randomisieren.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-adv-2
#| fig-align: center
#| fig-height: 5
#| fig-width: 7
#| fig-cap: "Randomized complete block design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` vollständig randomisiert in den vier Blöcken."
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = fert, cex = 1, shorten = "no", col = soil,
data = rcbd_fac3_book_tbl,
main = "Randomized complete block design (3-faktoriell)",
show.key = TRUE, key.cex = 1)
```
## Split plot design (3-faktoriell) {#sec-split}
Die Spaltanlage (eng. *split plot design*) ist eine häufig genutzte Variante, wenn wir eine Behandlungsfaktor nur in einem anderen Level eines zweiten Behandlungsfaktor randomisieren können. Wir können also eine Feldspur nur auf eine Art mechanisch bearbeiten und geben dann aber pro Feldspur verschiedene Dünger auf. Oder wir können in einen Stahl nur auf eine Art Futter zuführen, aber verschiedene Arten von Futter. Das gleiche kannst du dir auch mit einer Klimakammer vorstellen in der wir verschiedene Pflanzenlinien stellen können.
Beide Behandlungsfaktoren sind ineinander und dann natürlich im Block genestet. Hier aber erstmal das Modell mit den drei Faktoren für den besseren Überblick. Wir haben einmal die Düngung `fert`, dann den Faktor Boden `soil` sowie die verschiedenen Blöcke mit `block`.
$$
y \sim \overbrace{soil}^{f_1} + \underbrace{fert}_{f_2} + \overbrace{block}^{f_3}
$$
mit
- dem Faktor Dünger `fert` und den sechs Düngestufen als Level `fe_1` bis `fe_6`.
- dem Faktor Boden `soil` und den vier Bodenarten als Level `soil_1` bis `soil_4`.
- dem Faktor Block `block` und den zwei Leveln `block_1` und `block_2`.
In der @fig-mermaid-splitplot sehen wir die Abhängigkeitsstruktur des *split plot designs*. Wir haben den Faktor `fert` den wir in den Faktor `soil` nesten. Den Faktor `soil` nesten wir dann wiederrum in die Spalten `cols` des Blocks. Deshalb nennen wir das ja auch eine Spaltanlage. Ein Faktor ist immer in den Spalten angeordnet und der andere Faktor in der Spalte randomisiert. Wie immer wird es vielleicht klarer, wenn du dir dazu die @fig-exp-adv-split als Ergebnis des Versuchsdesigns anschaust.
```{mermaid}
%%| label: fig-mermaid-splitplot
%%| fig-width: 6
%%| fig-cap: "Abhängigkeitsstruktur des *split plot design*. Der Faktor `soil` ist in den Spalten `cols` der Blöcke randomisiert und der zweite Faktor `fert` innerhalb des anderen Faktors und somit auch in den Spalten."
flowchart LR
A(fert):::factor --- B(((nestet))) --> C(soil):::factor --- D(((nestet))) --> E(cols) --- F(block):::factor
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
```
::: callout-tip
## Modell zur Auswertung
Wir rechnen ein komplexeres gemischtes Modell mit allen möglichen zufälligen Effekten.
```{r}
#| eval: FALSE
fit <- lmer(yield ~ soil + fert + soil:fert +
(1|block) + (1|block:fert) + (1|block:soil),
data = split_tbl)
```
Oder können überlegen ein gemischtes Modell mit weniger zufälligen Effekten zu rechnen.
```{r}
#| eval: FALSE
fit <- lmer(yield ~ soil + fert + soil:fert +
(1|block),
data = split_tbl)
```
Die Entscheidung kannst du dann mit der [Modellselektion](https://jkruppa.github.io/stat-modeling-basic.html#sec-model-basic-compare) durchführen.
:::
Jetzt haben wir die volle Auswahl an Möglichkeiten. Ich zeige einmal wie man händisch das *split plot designs* erstellt. Dann schauen wir uns die Erstellung in `agricolae` einmal an und dann zeige ich noch die Variante in `dae`, die mich echt einiges an Zeit und Nerven gekostet hat. In sich ist das Paket `dae` ja logisch, aber die Dokumentation lässt für mich etwas zu wünschen übrig. Dennoch hier einmal die volle Dreifaltigkeit der Versuchsdesignerstellung, wenn ich mir schon die Mühe gemacht habe.
### ... mit `expand_grid()`
Ja, selber machen ist hier etwas mühsamer, aber wenn du die Schritte nachvollziehst, dann wird dir vermutlich das *split plot design* sehr viel klarer. Als erstes erschaffen wir die zwei Blöcke und die vier Bodenarten. Dabei sind in jedem Block die vier Bodenarten genestet. Dann gruppieren wir nach Block und randomiseren einmal die Bodenarten je Block. Im nächsten Schritt erweitern wir dann jede Block/Boden-Kombination um die sechs Düngerstufen. Dann gruppieren wir wieder, aber diesmal für jede Block/Boden-Kombination, um hier dann einmal die Düngestufen zu randomisieren. Dann lösen wir alle Gruppen auf und setzen fürunsere Faktoren dann noch die richtigen Labels.
```{r}
splitplot_long_tbl <- expand_grid(block = 1:2,
soil = 1:4) |>
group_by(block) |> # Gruppieren nach block
mutate(soil = sample(soil)) |> # Randomisieren von soil in block
expand_grid(fert = 1:6) |>
group_by(block, soil) |> # Gruppieren nach block und soil
mutate(fert = sample(fert)) |> # Randomisieren von fert in block und soil
ungroup() |>
mutate(fert = factor(fert, label = str_c("fe_", 1:6)),
block = factor(block, labels = str_c("Block ", 1:2)),
soil = factor(soil, labels = str_c("soil_", 1:4)))
splitplot_long_tbl
```
Da ich auch gleich wieder für die Blöcke separieren möchte, baue ich mir für jeden Block `rep` nochmal vier Spalten, je eine Spalte für eine Bodenbehandlung, und dann nochmal sechs Zeilen, eine für jede Düngestufe. Dann verbinden wir den randomisirten Datensatz mit den Pflanzgrid und schon können wir uns das *split plot design* einmal anschauen.
```{r}
splitplot_plot_tbl <- splitplot_long_tbl |>
bind_cols(expand_grid(rep = 1:2,
cols = 1:4,
rows = 1:6))
```
In der @fig-exp-adv-split siehst du einmal das *split plot design* dargestellt. Hier wird auch nochmal schön das Randomisierungsmuster klar. Du siehst auch warum das *split plot design* im Deutschen auch Spaltanlage heißt. Der Dünger ist in Spalten des Bodenfaktores randomisiert. Jeder Block hat also seine eigene Randomiserung der Böden.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-adv-split
#| fig-align: center
#| fig-height: 5
#| fig-width: 7
#| fig-cap: "Split plot design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` randomisiert in den zwei Blöcken. Der Faktor `fert` ist in dem Faktor `soil` genestet und randomisiert."
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = fert, cex = 1, shorten = "no", col = soil,
data = splitplot_plot_tbl ,
main = "Split plot design (3-faktoriell)",
show.key = TRUE)
```
### ... mit `{FielDHub}`
Für die Erstellung des *split plot design* in der Shiny App schaue dir die Anleitung unter [Generates a Split Plot Design mit `{FielDHub}`](https://didiermurillof.github.io/FielDHub/reference/split_plot.html) einmal an. Hier einmal der Code für die Erstellung in R. Wir haben jetzt vier Behandlungen für den Faktor `soil` sowie sechs Behandlungen für den Faktor `fe` vorliegen und wollen jeweils drei Blöcke anlegen. Jetzt einmal die Funktion `split_plot()` in der wir hier aber die Namen der Level übergeben können. Das ist dann wiederum eine Erleichterung zum dreifaktoriellen Design eben.
```{r}
split_plot_obj <- split_plot(wp = paste0("soil_", 1:4),
sp = paste0("fe_", 1:6),
reps = 3,
seed = 2023)
```
Dann können wir uns einmal das erstellte Feldbuch wiedergeben lassen, was du dann auch in Excel raus schreiben kannst. Wie immer brauchst du nicht alle Informationen. Auch hier würde ich das Feldbuch dann nutzen um mir selber in `{desplot}` die Abbildung zu erstellen. Die automatisierte Abbildungserstellung in `{FielDHub}` ist der manuellen weit unterlegen.
```{r}
split_plot_obj$fieldBook |>
as_tibble()
```
Dann erstellen wir auch schon die @fig-exp-basic-split-fieldhub mit dem *split plot design*. Leider ist die Abbildung nur noch okay. Die Einteilung mit dem `|` ist besser als beim faktoriellen Design, aber so richtig überzeugen tut mich die Darstellung nicht. Wir haben ja die wilden Farben, da wir jeweils die Faktorkombination wiedergeben. Leider helfen die Farben wenig und auch ist alles sehr unübersichtlich. Das Paket `{desplot}` kann definitiv mehr. Hier würde ich dann wirklich auf den Plot hier aus `{FielDHub}` verzichten und stattdessen mir den Plot wie oben in der @fig-exp-adv-split aus den den Feldbuch selber bauen.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-basic-split-fieldhub
#| fig-align: center
#| fig-height: 4
#| fig-width: 10
#| fig-cap: "Split plot design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` randomisiert in den zwei Blöcken. Der Faktor `fert` ist in dem Faktor `soil` genestet und randomisiert."
split_plot_obj |> plot()
```
### ... mit `{agricolae}`
Jetzt machen wir das Ganze nochmal mit dem R Paket `{agricolae}`. Hier müssen wir wieder die Vektoren mit den Faktoren und Leveln vorab erschaffen. Auch brauchen wir die Anzahl an Leveln für jeden Faktor. Am Ende randomisieren wir hier mal in vier Blöcke, einfach um es etwas anders zu machen.
```{r}
trt_fac1_soil <- str_c("soil_", 1:4)
n_trt_fac1_soil <- n_distinct(trt_fac1_soil)
trt_fac2_fert <- str_c("fe_", 1:6)
n_trt_fac2_fert <- n_distinct(trt_fac2_fert)
n_block <- 4
```
Für die Erstellung des *split plot design* nutzen wir die Funktion `design.split()`. Hier ist es wichtig, unbedingt die `serie = 0` als Option zu setzen. Wir brauchen die so entstehende Information in der Spalte `plots` sonst klappt es im Folgenden nicht die `plots` gleich den Spalten `cols` zu setzen. Ja, ich weiß, ist alles super suboptimal, aber so ist es eben in `{agricolae}`. Du musst noch einges machen, damit die Funktion auch einen visualisierbare Ausgabe wiedergibt.
```{r}
splitplot_obj <- design.split(trt1 = trt_fac1_soil,
trt2 = trt_fac2_fert,
r = n_block,
seed = 42, serie = 0)
```
Leider gibt es kein `sketch`-Objekt gerade hier wo man eins gebrauchen könnte. Deshalb also nochmal alles umbauen und die Zeilen `rows` und Spalten `cols` entsprechend ergänzen. Ich nenne jetzt nochmal ein wenig die Variablen um, einfach damit hier nicht immer das Gleiche passiert. Wichtig ist, dass wir hier immer für die Blöcke separieren. Wir müssen also in `desplot` definitiv die Option `| block` setzen, sonst klappt es mit der Darstellung nicht.
```{r}
splitsplot_book_tbl <- splitplot_obj$book |>
mutate(block = str_c("Block ", block),
cols = plots,
rows = as.numeric(splots)) |>
select(block, fert = trt_fac2_fert, soil = trt_fac1_soil, rows, cols)
```
In der @fig-exp-adv-3 siehst du das Ergebnis der Versuchsplanerstellung mit `{agricolae}`. Wir kriegen faktisch das Gleiche raus, gut die Faktoren sind anders permutiert, aber das wundert uns jetzt nicht. Am Ende musst du entscheiden, was dir besser gefällt. Je komplexer die Randomisierung ist, desto einfacher ist die Nutzung der `{agricolae}` Funktionen. Beim *split plot design* ist es so eine Sache, es lohnt sich bei sehr großen Anlagen dann schon.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-adv-3
#| fig-align: center
#| fig-height: 5
#| fig-width: 7
#| fig-cap: "Split plot design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` randomisiert in den vier Blöcken. Der Faktor `fert` ist in dem Faktor `soil` genestet und randomisiert."
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = fert, cex = 1, shorten = "no", col = soil,
data = splitsplot_book_tbl ,
main = "Split plot design (3-faktoriell)",
show.key = TRUE)
```
### ... mit `{dae}`
Auch das R Paket `{dae}` liefert die Möglichkeit ein *split plot design* zu erstellen. Die Idee von `dae` ist ja ein mehr generalisiertes Framework für die Erstellung von Versuchsplänen anzubieten. Deshalb gibt es ja nur die Funktion `designRandomize()`, die ein gegebenes Design nach feststehenden Regeln randomisiert. Die Idee ist gut, aber dafür musst du dich relativ tief in das Paket `dae` einarbeiten. Für mich war dieser kurze Code für das *split plot design* schon recht langwierig aus den Beispielen aus den [Notes on the use of dae for design](https://cran.r-project.org/web/packages/dae/vignettes/DesignNotes.pdf) herzuleiten. Prinzipiell ist es ja nicht schwer, aber es ist schon etwas mehr Auswand.
Wie funktioniert nun die Erstellung eines *split plot design*? Wir bauen uns erstmal unser Design der Positionen mit den Blöcken, Mainplots `MPlots` sowie den Subplot `SubPlots`. Dabei ist es wichtig zu wissen, dass die Subplots in den Mainplots genestet sind. Die Mainplots sind dann wiederum in den Blöcken genestet. Ja, hier hast du dann ein Paket spezifisches Naming der Optionen. Im zweiten Schritt packen wir dann den Faktor Boden `soil` und den Faktor `fert` zu den Positionen hinzu. Wichtig ist das `times = 4` am Ende, was nochmal die Faktorkombinationen von `soil` und `fert` jeweils in der Anzahl der Blöcke wiederholt. Wir haben ja hier vier Blöcke vorliegen, dass vergisst man ja mal schnell.
```{r}
split_sys <- cbind(fac.gen(list(Blocks = 4, MPlots = 4, SubPlots = 6)),
fac.gen(list(soil = str_c("soil_", 1:4),
fert = str_c("fe_", 1:6)), times = 4))
split_sys |>
as_tibble()
```
An der Ausgabe siehst du ganz schön, wie die Spalte Mainplots zu dem Faktor `soil` passt sowie die Spalte Subplots zu dem Faktor `fert`. Die Herausforderung ist nun, dass das R Paket `dae` davon ausgeht, dass du weißt, was ein *split plot design* ist und wie es aufgebaut ist. Es gibt also wenig Hilfestellung innerhalb der Funktionen. Dafür sind die Funktionen super flexibel und du kannst dir deine eigenen kreativen Designs bauen - wovon ich abrate.
Im nächsten Schritt nutzen wir die Funktion `designRandomize()` um die Faktoren `soil` und `fert` den Positionen in der Spalte `Block`, `MPlots` und `SubPlots` randomisiert zuzuweisen. Deshalb sind ja auch die Faktoren `soil` und `fert` der Option `allocated` (deu. *zugeteilt*) zugewiesen, da die beiden Faktoren ja den Spalten `Block`, `MPlots` und `SubPlots` als Empfänger (eng. `recipient`) zugeteilt werden. Soweit so gut, dass sagt aber noch nichts über das *Wie* der Zuteilung aus. Das machen wir dann mit der Option `nestet.recipients`. Der Option sagen wir jetzt wer in was genestet ist. Und durch diese Zuordnung erschaffen wir ein *split plot design*. Die Mainplots sind in den Blöcken genestet und dann die Subplots in den Mainplots und den Blöcken. Ja, da muss man erstmal drauf kommen. Hat bei mir gedauert, bis ich das hingekriegt habe.
```{r}
split_lay <- designRandomize(allocated = split_sys[c("soil", "fert")],
recipient = split_sys[c("Blocks", "MPlots", "SubPlots")],
nestet.recipients = list(MPlots = "Blocks",
SubPlots = c("MPlots", "Blocks")),
seed = 235805)
split_lay |>
as_tibble()
```
Und wir haben eine Randomisierung. Leider wissen wir nicht, ob die so geklappt hat, es ist ja wirklich schwer zusehen. Deshalb wollen wir uns die Ausgabe mal visualisieren. Das Gleiche Problem wie schon eben haben wir mit der Visualisierung mit `designGGPlot()`. Die Hilfeseite ist etwas unterkomplex für die Möglichkeiten der Funktion. So musste ich auch hier recht lange Rumspielen, bis ich die folgende @fig-exp-adv-4 erstellen konnte. Der Witz ist eigentlich hier, dass wir eine interne Funktion nutzen und auch das Objekt aus `designRandomize()` übergeben. Dennoch müssen wir uns selber mit den `row.factors()` und den `column.factors()` alles zusammenkleben. Ich habe es dann geschafft, aber es war mehr Zufall als Können. Das die Spalten die Mainplots sind, war mir klar, aber das die Zeilen durch die beiden Vektoren `Blocks` und `Subplots` so in der Reihenfolge gebildet werden war nur durch probieren zu lösen. Das eine aussagekräftige Legende fehlt, ist dann schon fast nicht mehr von Belang.
```{r}
#| echo: false
#| warning: false
#| label: fig-exp-adv-4
#| fig-align: center
#| fig-height: 6
#| fig-width: 6
#| fig-cap: "Split plot design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` randomisiert in den vier Blöcken. Der Faktor `fert` ist in dem Faktor `soil` genestet und randomisiert."
designGGPlot(split_lay, labels = "fert", cellalpha = 0.75, cellfillcolour.column = "soil",
row.factors = c("Blocks", "SubPlots"), column.factors = c("MPlots"),
blockdefinition = rbind(c(6, 1)),
blocklinecolour = "black", title = NULL)
```
Am Ende ermöglicht das Paket `{dae}` super flexibel experimentelle Design zu erstellen. Ich werde das Paket auch für Designs im Rahmen der Analyse von linearen gemischten Modellen nutzen. Dafür ist das Paket super. Du musst aber wissen, dass das Paket einiges an Einarbeitung benötigt und du auch die Philosophie hinter den Funktionen verstehen musst. Mal eben schnell, geht leider nicht. Dafür kann das Paket `{dae}` mehr als als andere Pakete im Bereich der Erstellung von experimentellen Designs.
## Subsampling {#sec-subsampling}
Die zentrale Idee des Subsamplings ist, dass wir nicht über die Pflanzen in unseren Blöcken mitteln. Wir nehmen also alle Pflanzen in unserem Block mit in die Analyse. Oder aber wir mitteln nicht das Gewicht der Ferkel für jede Bucht, sondern wollen die Gewichte individuell mit betrachten. Wir messen die Beobachtungen innerhalb der Behandlung und nehmen diese individuellen Beobachtungen auch mit ins Modell. Damit haben wir eigentlich das gleiche Modell wie auch bei dem *randomized complete block design* mit zwei Faktoren. Nur das wir eben in den Blöcken und Behandlungen nicht einen Mittelwert über alles haben, sondern eben sechs individuelle Werte für die sechs Beobachtungen in einer Behandlung/Block-Kombination in unserem Beispiel.
::: callout-tip
## Beispieldaten aus `{agridat}` zu Subsampling
Wenn dich nochmal Daten zum Subsampling interessiert dann schaue dir doch die Beispieldaten aus dem R Paket `{agridat}` an. Zum einen wäre der Datensatz [Rice RCB with subsamples](https://kwstat.github.io/agridat/reference/grover.rcb.subsample.html) sowie die Daten [Split-plot experiment of rice, with subsamples](https://kwstat.github.io/agridat/reference/gomez.splitplot.subsample.html) von Interesse. Du findest hinter dem Link dann auch die Hilfe wie du die Daten in R mit `data()` lädst.
:::
Das generelle Modell unterscheidet sich nicht von einem RCBD. Wir haben zwei Faktoren, einmal die Behandlung und einmal einen Block vorliegen. Damit sind wir noch sehr ähnlich. Den Unterschied macht gleich die Randomisierung und die Betrachtung der einzelnen individuellen Beobachtungen.
$$
y \sim \overbrace{trt}^{f_1} + \underbrace{block}_{f_2}
$$
mit
- dem Faktor Behandlung `trt` und den vier Behandlungsstufen als Level `ctrl`, `A`, `B` und `C`.
- dem Faktor Block `block` und den vier Leveln `block_1` bis `block_4`.
Es gibt jetzt zwei Arten, wie wir die Daten randomisieren und betrachten können. Zum einen randomisieren wir die Individuen vollständig auf die Behandlungen und die Behandlungen nesten wir dann in die Blöcke. Dann haben wir die Abhängigkeitssruktur in der @fig-mermaid-sub-1 vorliegen. Im Anschluss mitteln wir dann nicht über alle Beobachtungen. Prinzipiell ist es also ein sehr großes *randomized complete block design*.
```{mermaid}
%%| label: fig-mermaid-sub-1
%%| fig-width: 6
%%| fig-cap: "Subsampling vollständig innerhalb der Blöcke randomisiert. Die individuellen Beobachtungen sind über alle Behandlungen innerhalb der Blöcke randomisiert"
flowchart LR
A(trt):::factor --- B(((nestet))) --> C(block):::factor
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
```
Anders sieht es in der Variante zwei aus, hier haben wir die individuellen Beobachtungen in Clustern für eine einzelne Behandlung vorliegen. Die Beobachtungen sind innerhalb der Behandlungen genestet. In unserem Beispiel stehen die individuellen Pflanzen immer in einem Pflanztray. Du kannst dir aber auch Ferkel von verschiedenen Müttern in einer Bucht in einem Stall vorstellen. Es ergibt sich dann die @fig-mermaid-sub-2, wo unsere einzelnen Beobachtungen in den Behandlungen und dann in den Blöcken genestet sind.
```{mermaid}
%%| label: fig-mermaid-sub-2
%%| fig-width: 6
%%| fig-cap: "Subsampling vollständig innerhalb der Behandlung randomisiert. Die individuellen Pflanzen stehen in Trays und können nicht innerhalb der Behandlung `trt` randomisiert werden."
flowchart LR
A(tray) --- B(((nestet))) --> C(trt):::factor --- D(((nestet))) --> E(cols) --- F(block):::factor
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
```
Für die abschließende statistische Analyse kommt es darauf an, wie die individuellen Pflanzen in den Behandlungen stehen. Eventuell ist alles durcheinander gemischt, das heißt, die Pflanzen sind in den Behandlungen und Blöcken randomisiert. Dann ist es recht einfach mit einem linearen Modell zu machen.
Wenn du aber die Pflanzen in einem Pflanztray oder irgendwie anders fixierst hast, dann macht es schon einen Unterschied. In dem Fall sind die Pflanzen ja individuell nicht unabhängig voneinander. Die Pflanzen sind ja in einem Tray zusammen. Dann nutzen wir ein lineares gemischtes Modell.
Im Folgenden generieren wir einmal beide Fälle, wobei der zweite Fall - Pflanzen stehen zusammen - als Subsampling angesehen wird. Häufig kann dies auch bei Tierversuchen passieren, wenn zum Beispiel Ferkel in einer Bucht in einem Stahl stehen. Da aber ab und zu beides vorkommt, hier dann einmal beides. Wir können das Subsampling nur selber mit `expand_grid()` durchführen, da weder das R Paket `{agricolae}` noch das R Paket `{dae}` eine entsprechende Implementierung der Funktionalität zur Generierung eines Datensatzes mit Subsampling aufweist.
::: callout-tip
## Modell zur Auswertung
Vollständig im Block randomisiert, wir rechnen ein lineares Modell.
```{r}
#| eval: FALSE
fit <- lm(drymatter ~ trt + block + trt:block, data = subsampling_tbl)
```
Behandlung im Block randomisiert, dass heißt, dass die Pflanzen zusammen in einem Tray stehen, wir rechnen ein gemischtes Modell. Der Tray-Effekt wird dann im zufälligen Effekt `(1 | trt:block)` modelliert.
```{r}
#| eval: FALSE
fit <- lmer(drymatter ~ trt + block + (1 | trt:block), data = subsampling_tbl)
```
:::
Das Subsampling geht nur im Selbermachen, also machen wir das dann einmal mit der Funktion `expand_grid()`. Als erstes erstellen wir einmal ein Design mit vier Blöcken und vier Behandlungen. Dann ergänzen wir noch pro Behandlung sechs individuelle Beobachtungen. Am Ende randomisieren wir einmal über die Behandlungen in den Blöcken. Wichtig ist hier, dass wir die Daten nach den Blöcken gruppieren müssen und schon vorher alle individuellen Pflanzen generieren.
```{r}
subsampling_long_tbl <- expand_grid(block = 1:4,
trt = 1:4,
rep = 1:6) |>
mutate(trt = factor(trt, labels = c("ctrl", "A", "B", "C")),
block = factor(block, labels = str_c("Block ", 1:4)),
pid = 1:n(),
pid_text = str_pad(1:96, width = 2, pad = "0")) |>
group_by(block) |>
mutate(trt = sample(trt)) ## Randomisierung in den Blöcken
```
Dann bauen wir noch unser Grid auf. Wir wollen insgesamt sechzehn Spalten, dass sind die Blöcke/Behandlungs-Kombinationen. Dann brauchen wir noch jeweils sechs Zeilen für die individuellen Pflanzen. Dann können wir die Daten schon nehmen und in der @fig-exp-adv-10-1 visualisieren. Bevor wir die Abbildung besprechen, generieren wir uns nochmal die anderen Daten, wo die Behandlungen in den Blöcken randomisiert werden.
```{r}
subsampling_plot_tbl <- subsampling_long_tbl |>
bind_cols(expand_grid(cols = 1:16, rows = 1:6))
```
Eine andere Möglichkeit ist, dass unsere Pflanzen in einem Pflanztray stehen. Damit haben wir keine Randomisierung mehr in den Blöcken, sondern wir können nur unsere Behandlungen in den Blöcken randomisieren. Die Positionen in den Trays sind ja fix. Daher Randomisieren wir erst den Faktor der Behandlungen in den gruppierten Blöcken. Dann ergänzen wir noch die sechs individuellen Pflanzen.
```{r}
subsampling_tray_long_tbl <- expand_grid(block = 1:4,
trt = 1:4) |>
mutate(trt = factor(trt, labels = trt_fac1_soil),
block = factor(block, labels = str_c("Block ", 1:4))) |>
group_by(block) |>
mutate(trt = sample(trt)) |> # Randomisierung der Behandlungen in den Blöcken
expand_grid(rep = str_pad(1:6, width = 2, pad = "0"))
```
Jetzt bauen wir uns noch das Grid für die Positionen zusammen, da bleibt alles beim alten. Wir brauchen auch wieder sechzehn Spalten, dass sind die Blöcke/Behandlungs-Kombinationen. Dann noch die sechs Zeilen für das Pflanzentray. In der @fig-exp-adv-10-2 sehen wir einmal das Ergebnis der Randomisierung. Alle Pflanzen sind jetzt in einer Spalte randomisiert, die eben das Tray der Pflanzung entspricht.
```{r}
subsampling_tray_plot_tbl <- subsampling_tray_long_tbl |>
bind_cols(expand_grid(cols = 1:16, rows = 1:6))
```
In der @fig-exp-adv-10 sehen wir nochmal die beiden Arten des Subsamplings miteinander vergleichen. In der @fig-exp-adv-10-1 sind alle Behandlungen mit ihren sechs Wiederholungen zufällig im Block verteilt. Im anderen Fall sehen wir in der @fig-exp-adv-10-2, wenn die Pfanzen zusammen in einem Tray stehen und alle die gleiche Behandlung kriegen. Es können aber auch Schweine in einer Bucht sein und die Blöcke sind dann die einzelnen Ställe.
```{r}
#| echo: false
#| warning: false
#| label: fig-exp-adv-10
#| fig-align: center
#| fig-height: 5
#| fig-width: 7
#| fig-cap: "Verschiedene Arten des Subsamplings. Einmal vollständig innerhalb der Blöcke randomsiert und einmal die Behandlung randomisiert."
#| fig-subcap:
#| - "... vollständig innerhalb Block randomisiert (`subsampling_plot_tbl`)."
#| - "... Behandlungen innerhalb Block randomisiert (`subsampling_tray_plot_tbl`)."
#| layout-nrow: 1
#| column: page
desplot(data = subsampling_plot_tbl, flip = TRUE,
form = block ~ cols + rows | block,
out1 = trt, out1.gpar = list(col = "grey", lty = 1),
out2 = pid, out2.gpar = list(col = "grey", lty = 1),
main = "Subsampling (2-faktoriell, vollständig randomisiert)",
text = pid_text,
cex = 1, show.key = FALSE, shorten = "no", col = trt)
desplot(data = subsampling_tray_plot_tbl, flip = TRUE,
form = block ~ cols + rows | block,
out1 = trt, out1.gpar = list(col = "grey", lty = 1),
#out2 = pid, out2.gpar = list(col = "grey", lty = 1),
main = "Subsampling (2-faktoriell, Behandlung randomisiert)",
text = rep,
cex = 1, show.key = FALSE, shorten = "no", col = trt)
```
## Incomplete block design {#sec-incomplete}
Bei dem *incomplete block design* haben wir eigentlich ein normales Design, wie das 3-faktorielle *randomized complete block design* vorliegen, aber die Randomiserung ist nicht vollständig über alle Blöcke. Wir nehmen einmal die beiden Spezialfälle [Alpha design](#sec-alpha) und [Augmented design](#sec-augment) aus der Betrachtung raus. Wir schauen uns hier nochmal zu Vervollständigung das *randomized complete block design* mit nicht kompletten Blöcken an. Ich persönlich stehe dem Design etwas zwiegesaplten gegenüber. Ich würde lieber die Anzahl der Gruppen der Behandlungen reduzieren als die Blöcke. Wenn du ein Pilotversuch machst um zu schauen welche Linien wachsen oder aber ob es überhaupt einen Effekt einer Behandlung gibt, dann mag das *incomplete block design* einen Sinn machen. Dann sparst du wirklich Ressourcen und Zeit.
::: callout-caution
## Verliere keine Behandlung im *incomplete block design*
Achte immer darauf, dass du genug Blöcke hast um auch alle Behandlungslevel aufzunehmen. Schaue lieber doppelt, ob du auch wirklich alle Behandlungen in deinem Versuchsplan drin hast. Es kann schnell passieren, dass du auf einmal nur ein Level der Behandlung oder gar keins mehr in deinem Design vorliegen hast. Das R Paket `agricolae` versucht das *incomplete block design* so zu bauen, dass du keine Behandlungen verlierst.
:::
Wir immer erstmal unser Modell mit den zwei Behandlungsfaktoren `soil` und `fert`. Dann kommt noch der nicht komplette Block `inc_block` hinzu.
$$
y \sim \overbrace{soil}^{f_1} + \underbrace{fert}_{f_2} + \overbrace{inc\_block}^{f_3}
$$
mit
- dem Faktor Dünger `fert` und den sechs Düngestufen als Level `fe_1` bis `fe_6`.
- dem Faktor Boden `soil` und den vier Bodenarten als Level `soil_1` bis `soil_4`.
- dem Faktor Block `block` und den zwei Leveln `block_1` und `block_2`.
In der @fig-mermaid-incomplete-1 sehen wir nochmal die Abhängigkeitsstruktur des *incomplete block design*. Die beiden Faktoren `soil` und `fert` sind in den Blöcken genestet. Wir stellen hier die Linien als gestrichelt dar um auszudrücken, dass wir nur eine nicht komplette Zuordnung der Behandlungsfaktoren zu den Blöcken haben.
```{mermaid}
%%| label: fig-mermaid-incomplete-1
%%| fig-width: 6
%%| fig-cap: "Schema der Abhängigkeitsstruktur des *incomplete block design* für das 3-faktorielle *randomized complete block design*. Die gestrichelte Linie stellt die nicht vollständige Zuordnung der Behandlungsfaktoren zu den Blöcken dar. Daher sind die Blöcke `inc_block` nicht komplett."
flowchart LR
A(fert):::factor -.- B(((nestet))) -.-> E(inc_block):::inc_block
C(soil):::factor -.- D(((nestet))) -.-> E
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
classDef inc_block fill:#CC79A7,stroke:#333,stroke-width:0.75px
```
Dann das Ganze nochmal in dem R Paket `{agricolae}` mit einem 2-faktoriellen *randomized complete block design* mit dem Fall der nicht kompletten Blöcke. Wir können in der Standardfunktion von `agricolae` nur ein 2-faktorielles *incomplete block design* abbilden. Deshalb dieses etwas kleinere Modell mit einem Behandlungsfaktor `trt` und einem Block `block`. Wir werden dann die Behandlungsgruppen nicht komplett auf alle Blöcke randomisieren.
$$
y \sim \overbrace{trt}^{f_1} + \underbrace{inc\_block}_{f_2}
$$ mit
- dem Faktor Behandlung `trt` und den vier Behandlungsstufen als Level `ctrl`, `A`, `B` und `C`.
- dem Faktor Block `soil` und den drei Blöcken mit den Leveln `1`, `2`, `3`.
Dann schauen wir uns nochmal die die Abhängigkeitsstruktur des *incomplete block design* in der @fig-mermaid-incomplete-2 an. Wir sehen, dass Design ist relativ einfach. Die Behandlungen sind in den Blöcken genestet, aber wir haben nicht alle Behandlungen in allen Blöcken vorliegen.
```{mermaid}
%%| label: fig-mermaid-incomplete-2
%%| fig-width: 6
%%| fig-cap: "Schema der Abhängigkeitsstruktur des *incomplete block design* für das 2-faktorielle *randomized complete block design*. Die gestrichelte Linie stellt die nicht vollständige Zuordnung der Behandlungsfaktoren zu den Blöcken dar. Daher sind die Blöcke `inc_block` nicht komplett."
flowchart LR
C(trt):::factor -.- D(((nestet))) -.-> E(inc_block):::inc_block
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
classDef inc_block fill:#CC79A7,stroke:#333,stroke-width:0.75px
```
Wir oben schon geschrieben, die Verwendung von dem R Paket `agricolae` schützt dich davor, keine der Behandlungslevel zu verlieren. Bei dem Selbermachen musst du dann immer schauen, ob dir der Versuchsplan so passt und alle Behandlungen einigermaßen vertreten und randomisiert sind.
::: callout-tip
## Modell zur Auswertung
Wir rechnen ein multiples lineares Modell für die statistische Analyse. Vermutlich reicht die Fallzahl nicht für alle Interaktionen zwischen allen Faktoren.
```{r}
#| eval: FALSE
fit <- lm(drymatter ~ fert + soil + fert:soil + inc_block, data = inc_block_tbl)
```
Wir rechnen ein multiples lineares Modell für die statistische Analyse.
```{r}
#| eval: FALSE
fit <- lm(drymatter ~ trt + inc_block + trt:inc_block, data = inc_block_tbl)
```
:::
### ... mit `expand_grid()` (3-faktoriell)
Bei dem Selbermachen ist es jetzt eigentlich wie bei dem 3-faktoriellen *randomized complete block design*. Deshalb schaue einfach nochmal oben nach, wenn dir Teile der Generierung unklar sind. Hier und dort ändert sich ja bei mir immer nur leicht der Code, damit du auch mal was anderes siehst. Es gibt ja viele Wege nach Rom.
```{r}
incomplete_long_tbl <- expand_grid(block = 1:2, soil = 1:4, fert = 1:6) |>
mutate(block = factor(block, labels = str_c("Block ", 1:2)),
soil = factor(soil, labels = str_c("soil_", 1:4)),
fert = factor(fert, label = str_c("fe_", 1:6))) |>
group_by(block) |>
mutate(fert = sample(fert),
soil = sample(soil))
incomplete_long_tbl
```
Noch haben wir kein *incomplete block design* vorliegen, dass bauen wir uns jetzt. Wir rasieren nämlich einfach die letzten beiden Zeilen unseres *randomized complete block design* ab. Damit haben wir dann nur noch ein $4\times4$-Block übrig. Sonst wäre es ja ein Block mit sechs Zeilen.
```{r}
incomplete_plot_tbl <- incomplete_long_tbl |>
bind_cols(expand_grid(rep = 1:2, cols = 1:4, rows = 1:6)) |>
filter(rows <= 4)
```
Jetzt müssen wir uns die Sachlage einmal anschauen. In der @fig-exp-adv-7 sehen wir einmal den Versuchsplan für das *incomplete block design*. Du musst immer schauen, ob wir wirklich alle Level jedes Behandlungsfaktors auch vorliegen haben und ob sich die Level auch einigermaßen gut verteilen. Das *incomplete block design* ist da etwas frickelig, ob es dann auch wirklich gut ist. Aber du kannst auch nciht alle Kombinationen vorliegen haben, denn sonst wäre es ja auch kein *incomplete block design*.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-adv-7
#| fig-align: center
#| fig-height: 5
#| fig-width: 7
#| fig-cap: "Incomplete block design (3-faktoriell) für die Faktoren Dünger `fert` sowie `soil` randomisiert in den zwei Blöcken."
desplot(block ~ cols + rows | block, flip = TRUE,
out1 = rows, out1.gpar = list(col = "grey", lty = 1),
out2 = cols, out2.gpar = list(col = "grey", lty = 1),
text = fert, cex = 1, shorten = "no", col = soil,
data = incomplete_plot_tbl,
main = "Incomplete block design (3-faktoriell)",
show.key = TRUE)
```
### ... mit `{FielDHub}`
Für die Erstellung des *incomplete block design* in der Shiny App schaue dir die Anleitung unter [Generates a Resolvable Incomplete Block Design mit `{FielDHub}`](https://didiermurillof.github.io/FielDHub/reference/incomplete_blocks.html) einmal an. Hier einmal der Code für die Erstellung in R. Leider ist das ganz viel rumprobieren, bis es dann mit den Parametern klappt. Dehalb empfhele ich dir hier auf jeden Fall die Shiny App, da kannst du viel einfach rumprobieren und mit den Anzahlen an Behandlungen und Blöcken spielen. Hier generierst du sehr viele Fehlermeldungen und versuchst dann ein Design zu finden was passt, dann musst du es ja noch visualisieren. Die Schritte sind dann in der Shiny App besser miteinander verwoben.
```{r}
incomplete_block_obj <- incomplete_blocks(t = 6,
k = 3,
r = 3,
seed = 1984)
```
Dann können wir uns einmal das erstellte Feldbuch wiedergeben lassen, was du dann auch in Excel raus schreiben kannst.
```{r}
incomplete_block_obj$fieldBook |>
as_tibble()
```
Dann erstellen wir auch schon die @fig-exp-basic-incomplete-fieldhub mit dem *split plot design*. Da wir hier wieder eine relativ einfache Abbildung haben, geht es dann auch mit der Darstellung in `{FielDHub}`. Leider sind auch hier die Zahlen teilweise sehr schwer auf den Farben zu erkennen.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-basic-incomplete-fieldhub
#| fig-align: center
#| fig-height: 4
#| fig-width: 5
#| fig-cap: "Incomplete block design (2-faktoriell) für den Faktor `trt` randomisiert in den sechs Blöcken."
incomplete_block_obj |> plot()
```
### ... mit `{agricolae}` (2-faktoriell)
Der Vorteil von `{agricolae}` und der Funktion `design.bib()` ist ganz klar die Effizienz und die optimale Verteilung der Behandlungsfaktoren auf die Blöcke. Wir können hier sicher sein das optimale Ergebnis der Randomisierung zu erhalten und müssen auch keine Sorge haben, das Level der Behandlungen abhanden gekommen sind. Deshalb kann ich hier die Funktion `design.bib()` in `{agricolae}` sehr empfehlen.
Die Namen der Optionen sind nicht die Stärke von `{agricolae}`, also haben wir hier wieder die alten Probleme. Wir müssen darauf achten, dass wir die Anzahl an Blöcken mit `k` angeben. Wenn du zu wenig Blöcke bei zu vielen Behandlungslevel wählst, dann gibt dir die Funktion einen Fehler wieder und kein Ergebnis. An der Ausgabe der Funktion kannst du dann gleich deine Effizienz gegenüber einer vollständigen Randomisierung abschätzen.
```{r}
#| message: false
#| warning: false
bib_obj <- design.bib(trt = c("ctrl", "A", "B", "C"),
k = 3, seed = 543, serie = 2)
```
Wichtig ist hier, dass wir auf den `Efficiency factor` schauen. Wir sind mit unserer Randomisierung des *incomplete block design* ungefähr zu 89% so gut wie das gleichwertige *randomized complete block design*. Das ist schonmal gut. Wir sparen Platz auf dem Feld oder dem Stall und sind trotzdem nicht so viel schlechter als ein vollständig randomisiertes Design. Aber Achtung, die Aussagekraft ist hier geringer, für Pilotstudien ist das immer okay, für echte Feldstudien zur Zulassung oder für wissenschaftliche Publikationen eher weniger bis gar nicht.
Wir müssen in diesem Fall nicht so viel an den Daten bereinigen um uns gleich die Visualisierung anzuschauen. Das geht hier relativ einfach. Wir benennen nur die Spalten um und sind schon soweit fertig.
```{r}
bib_book_tbl <- bib_obj |>
pluck("book") |>
as_tibble() |>
set_names(c("plots", "block", "trt"))
bib_book_tbl
```
Dann ergänzen wir noch die Informationen zu den Spalten und den Zeilen, so dass wir dann auch gleich `desplot()` nutzen können. Wir müssen hier nochmal nach den Spalten gruppieren, damit wir die Zeilen sauber zugeordnet kriegen. Wir haben ja nicht für jede Behandlung eine Zeile sondern eben weniger. Ja, man hätte auch mit der Information von `k = 3` hier arbeiten können. So finde ich es aber schöner.
```{r}
bib_plot_tbl <- bib_book_tbl |>
mutate(cols = as.numeric(block)) |>
group_by(cols) |>
mutate(rows = 1:n())
```
Einmal der Vergleich zum `sketch` aus der Ausgabe von `design.bib()`. Siehst gut aus und entspricht unseren Erwartungen. Du siehst hier gut, dass im ersten Block die Kontrolle gar nicht enthalten ist.
```{r}
bib_obj |>
pluck("sketch") |>
t()
```
In der @fig-exp-adv-8 sehen wir die Visualisierung unseres *incomplete block design* für den Faktor `trt` randomisiert in den vier Blöcken. Ich finde die farbliche Hinterlegung sehr schön, da sieht man viel schneller die eigentliche Zuordnung der Behandlungslevel zu den Blöcken. Auch siehst du schön, dass in dem vierten Block die Behandlung `C` gar nicht vorkommt. Was du dann nicht misst, kann auch nicht in der statistischen Auswertung berücksichtigt werden.
```{r}
#| message: false
#| warning: false
#| label: fig-exp-adv-8
#| fig-align: center
#| fig-height: 3
#| fig-width: 6
#| fig-cap: "Incomplete block design (2-faktoriell) für den Faktor `trt` randomisiert in den vier Blöcken."
desplot(trt ~ cols + rows, flip = TRUE,
text = trt, cex = 1, shorten = "no",
out1 = block,
data = bib_plot_tbl,
main = "Incomplete block design (2-faktoriell)",
show.key = FALSE)
```
## Strip plot design (3-faktoriell) {#sec-strip-plot}
Das Streifenparzellenversuch (eng. *strip plot design*) ist eigentlich ein doppeltes *split plot design* oder eben eine Erweiterung. Wo das *split plot design* in Spalten orientiert ist ergänzen wir jetzt durch das *strip plot design* noch die horizontale Dimension. Das heißt wir haben eine Randomisierung und Abhängigkeitsstruktur in den Spalten sowie in den Zeilen. Dabei orientieren wir uns an dem [Strip-plot experiment of rice](https://kwstat.github.io/agridat/reference/gomez.stripplot.html) aus @gomez1984statistical. Du findest dann das Datenbeispiel auch in `{agridat}` als `data(gomez.stripplot)`. Ich habe versucht jetzt die Daten in der Form nachzubauen und somit kannst du dann auch ähnliche Experimente planen.
Unser *strip plot design* hat drei Wiederholungen mit `block`, eine Sorte `trt_gen` als horizontaler Streifen und Stickstoffdünger `trt_nitro` als vertikaler Streifen. Das ganze klingt immer super abstrakt, deshalb gehe gerne gleich zur Auflösung in die @fig-exp-adv-9. Manchmal ist es besser sich gleich die Visualisierung anzuschauen und dann das folgende Modell.
$$
y \sim \overbrace{trt_{gen}}^{f_1} + \underbrace{trt_{nitro}}_{f_2} + \overbrace{block}^{f_3}
$$
mit
- dem Faktor Behandlung `trt_gen` und den sechs Linien als Level `G_1` bis `G_6`.
- dem Faktor Dünger `trt_nitro` und den drei Düngestufen als Level `0`, `60` und `120`.
- dem Faktor Block `block` und den drei Leveln `1` bis `3`.
In der @fig-mermaid-stripplot dann einmal die Abhängigkeitsstruktur des *split plot design*. Wichtig ist hierbei, dass die genetischen Linien `trt_gen` in den Zeilen `rows` der Blöcke genestet und randomisiert sind. Wir sehen aber die Zeilen später nicht in der Auswertung. Die Zeilen werden durch den Faktor Block repräsentiert. Genauso ist es mit der Stickstoffbehandlung `trt_nitro`, die dann in den Spalten `cols` des Blocks genested und randomisiert wird. Auch hier sehen wir dann in der späteren Auswertung nur den Block.
```{mermaid}
%%| label: fig-mermaid-stripplot
%%| fig-width: 6
%%| fig-cap: "Schema der Abhängigkeitsstruktur des *split plot design*. Die genetischen Linien `trt_gen` sind in den Zeilen `rows` der Blöcke genestet und randomisiert. Die Stickstoffbehandlung `trt_nitro` ist dann in den Spalten `cols` des Blocks genested und randomisiert."
flowchart LR
A(trt_gen):::factor --- B(((nestet))) ---> C(rows) --- E(block):::factor
F(trt_nitro):::factor --- G(((nestet))) ---> H(cols) --- E
classDef factor fill:#56B4E9,stroke:#333,stroke-width:0.75px
```
::: callout-tip
## Modell zur Auswertung
Wir rechnen hier ein lineares gemischtes Modell für die Auswertung der Daten.
```{r}
#| eval: FALSE
fit <- lmer(yield ~ trt_gen + nitro + trt_gen:nitro +
(1 | block) + (1 | block:nitro) +
(1 | block:trt_gen),
data = stripplot_tbl)
```
:::
### ... mit `{agricolae}`
Das *strip plot design* kann ich selber machen, aber hier fängt es dann an wirklich so kompliziert zu werden, dass es sich nicht mehr lohnt. Dafür haben wir dann in `{agricolae}` die passende Funktion. Die würde ich dann auch immer empfehlen zu nutzen, dass macht das Leben wirklich einfacher. Wir müssen zwar wieder erstmal unsere Faktoren als Vektoren bauen und dann alles in den passenden Objekten ablegen, aber das ist nur wenig Aufwand im Vergleich zum Selbermachen.
```{r}
trt_gen <- str_c("G", 1:6)
n_trt_gen <- n_distinct(trt_gen)
trt_nitro <- c(0, 60, 120)
n_trt_nitro <- n_distinct(trt_nitro)
n_block <- 3
```
Jetzt kommt der einzige wichtige Teil, damit du später die Zeilen- und Spaltenpositionen einfach zuweisen kannst. Du musst in die Option `trt1 =` den *kürzen* Vektor packen. Hier hat der Dünger nur drei Level, also ist der Dünger in der Option `trt1 =` vertreten. Dann nehmen wir den längeren Vektor für die genetischen Linien in die zweite Option `trt2 =` mit rein. Nimm also erst den kurzen Vektor dann den langen Vektor, dann klappt es auch mit der späteren Zuweisung der Positionen leichter.
```{r}
strip_obj <- design.strip(trt1 = trt_nitro, # kürzerer Vektor (hier der Länge 3)
trt2 = trt_gen, # längerer Vektor (hier der Länge 6)
r = n_block, serie = 2, seed = 543)
```
Für die Zuweisung der richtigen Zeilen `rows` und Spalten `cols` im Positionsgrid müssen wir dann nicht so viel machen, wenn wir oben darauf geachtet haben, die Behandlungsfaktoren in der richtigen Reihenfolge den Optionen zuzuweisen. Die Zuweisung der `cols` funktioniert also nur, wenn du in der Funktion `design.strip()` wirklich den kürzeren vor dem längeren Vektor zuweist.
```{r}
strip_plot_tbl <- strip_obj$book |>
bind_cols(expand_grid(cols = 1:(n_trt_nitro*n_block),
rows = 1:n_trt_gen))
```
Dann sind wir schon soweit uns einmal in der @fig-exp-adv-9 das *strip plot design* anzuschauen. Du kannst sehr schon sehen, wie die genetischen Linien in den Zeilen randomisiert sind und die verschiedenen Stickstoffgaben dann in den Zeilen. Die drei Blöcke werden dann durch die dicken Linien dargestellt. Da wir es bei dem *strip plot design* häufig mit Großversuchen zu tun haben sind die Blöcke meistens durch Feldwege oder andere bauliche Bedingungen voneinander getrennt.
```{r}
#| message: false