Several factors are believed to influence the rate and degree of “great success” when producing items using crafting skills in UWO. However, conflicting opinions exist as to what some of the different factors exactly do. Here I experimentally investigate the effects of three factors on a) the probability of great success (i.e. how often), and b) the degree of great success (i.e. how much improvement). The study was performed making the weapon Gorz – a rank 17 casting recipe where a great success can result in both increased attack and/or defence, as well as durability. In total 2786 Gorz items were made under controlled conditions, 1295 of which were great success. I conclude that the special production EX and the item improvement Oxford skill influence neither the rate nor the degree of great success. Instead, these features may exclusively apply to recipes that improve an already existing item (known as forging). In contrast, having the casting skill refined hugely increased the degree of great success, but had no influence on the rate of great success. Without refined skill the maximum improvement value was +40% of the base value, after refining the skill this increased to +95% for attack and defence.

Keywords: Uncharted Waters Online; Great Success; Casting; UWO myths.

1 Introduction

Most recipes for crafting in UWO have at least two possible outcomes: normal success and great success. Some also have a failure possibility resulting in nothing being produced (or a qualitatively inferior item). Some recipes result in a superior item with another name when obtaining a great success, other recipes keep the same name. (In general, name-changed items are qualitatively different ones, such as for example master’s cannons). In both cases, some of the stats of items from such great success, such as durability, attack or defence value have been improved compared to normal success. Many recipes have different degrees of great success, some in clearly defined tiers such as penetration power on cannons. The goal of crafting is often to produce an item with the highest possible tier of some stat(s), in as few attempts as possible. It is therefore of much interest to understand how to increase both the overall rate of great success and how to increase the degree of great success.

From previous more or less anecdotal experience we know with some level of certainty about several factors that influence the rate of success. Apart from using the Master’s Secrets item (which will give 100% great success), this includes:

Production success 1–4 from Oxford
Producing at sea with workshop ship skill
Having the relevant Meister title on (requires r15 in the production skill) with apprentice nearby
Having the Magician tarot card active

Other factors believed to influence the rate and/or the degree of great success include:

Skill rank
Refined skill
Item improvement tech from Oxford
Sailor Equipment for workshop
Special production EX
Aide in paymaster job
Paymaster aide trait value

Confusion and uncertainty exist about what the exact effect of some of these factors is. It has been unclear what influences the rate of great success and what influences the degree.

Some of these factors were studied in a previous report (Stinker 2023) in the context of the outcome of advanced alchemy experiments, which may or may not be directly applicable to normal production recipes. In that study I found that skill rank influences the degree of success, but not the rate, and no effect of Oxford skills, paymaster aide or EX equipment. However, this is complicated by the fact that those alchemy experiments have both failure, normal success and “huge” success, and that the normal success comes in qualitatively different tiers – it was only the tiers of normal success that was influenced by skill rank, not the probability of “huge” success.

In order to test the effect of some of these factors on success of ordinary crafting recipes, I conducted controlled experiments by manipulating certain factors while keeping everything else constant. In particular, I wanted to test the effect of the special production EX and the item oxford tech skill, since conflicting claims exist about what these do. I was also interested in how large is the effect of having refined the relevant crafting skill. I tested the effect of the following factors on both the rate and degree of great success:

Having refined the casting skill or not
Oxford Item Improvement active or not
Special production EX equipped or not (r6, Golden Ark)

I focused on the r17 casting recipe of making Gorz, fixed at the craftsman in Calcutta (Figure 1). Gorz is a mace weapon (originally known from Persian mythology, later introduced to India) that as base has 50 attack, 10 defense, 50 durability, and a great success is signified in-game by a statement that you have produced a superior item. Earlier experience revealed that all three stats can be increased by great success. Thus, several parameters can be measured at the same time, giving more data and also an opportunity to study if these are related. It was therefore deemed an ideal recipe, and the needed materials (steel and iron) are also easy to make from raw materials sold nearby.

FIGURE 1: The recipe for making Gorz at the Craftsman in Calcutta

One disadvantage of the Gorz recipe is that the required casting rank is high, so no effort was made to test the effect of skill rank. The effect of skill rank should instead be tested with a low rank recipe, where one may expect a more pronounced (and therefore easier to measure) effect of skill rank difference between being only at the required rank vs. at the maximum rank.

2 Methods

The stats of the weapons was recorded as each were made (normal or great success; durability; attack power; defence power). The weapons were made in batches of 20 without closing the production window. After a batch was done, entered values were double checked and then the produced weapons were disposed of to make space for a new batch. In a few cases less than 20 were produced in a batch because I ran out of materials. These incomplete batches were included in the analyses of degree of great success, and some in analyses about rate of success, but excluded from analyses of rate of success on batch level.

For each treatment group I aimed to make at least 250 great success weapons, instead of trying to keep the total sample size the same. Thus, if there was any difference in the rate of success the totals might differ as a consequence.

Due to fatigue I did not make all 8 possible treatments (combinations of the 3 factors with 2 levels each, \(2^3\) = 8), but if we ignore interaction effects this is also not necessary. Data was obtained for the 5 treatments shown in Table 1.

Treatment	Factors			Batches		Total Produced
TABLE 1: Experimental treatment groups and sample sizes
Treatment	Casting	Oxford	Ex	Complete	Incomplete	Total Produced
A	Refined	Yes	Yes	26	1	530
B	Refined	Yes	No	28	1	576
C	Refined	No	No	29	-	580
D	Unrefined	Yes	Yes	28	-	560
E	Unrefined	Yes	No	27	-	540
F (not done)	Unrefined	No	Yes	-	-	-
G (not done)	Unrefined	No	No	-	-	-
H (not done)	Refined	No	Yes	-	-	-
Total	—	—	—	138	2	2786

All experiments were done with:

Aide in paymaster job with 120SS paymaster traits
Oxford production success 1,2,3 and 4 active
Rank 20 casting
Casting Meister title active and apprentice nearby
No tarot active

These settings were chosen as to potentially give the highest possible rate of great success while at the same time keeping conditions similar in all respects apart from the manipulated factors. However, the experiments with refined skill were done with another character than the ones with unrefined skill (level 95/95/95 vs 51/77/44, respectively). The experiments with refined casting were done in October 2024 (Shipmate patch), and unrefined in November 2024 (Summer Triangle patch).

This report was written using R Markdown in RStudio. R code used to generate the results can be found in the Appendix together with software versions used.

3 Results

3.1 Rate of success

3.1.1 Overall and per treatment

In total, I made 2786 Gorz of which 1295 were a great success. Thus the overall rate of great success was 46.5%. Table 2 shows the proportions of normal and great success for Refined/Unrefined (\(\chi^{2}\) = 0.03, df = 1, p = 0.86), Oxford/No Oxford (\(\chi^{2}\) = 0.72, df = 1, p = 0.39) and EX rank 6/No EX (\(\chi^{2}\) = 1.16, df = 1, p = 0.28).

Comparison	Success				Total
TABLE 2: Proportion great success. Note that the the totals are the same in all the 3 comparisons
Comparison	Great		Normal		Total
Casting
Refined	781	(46.3%)	905	(53.7%)	1686
Unrefined	514	(46.7%)	586	(53.3%)	1100
Total	1,295	—	1,491	—	2,786
Oxford
Oxford	1035	(46.9%)	1171	(53.1%)	2206
No Oxford	260	(44.8%)	320	(55.2%)	580
Total	1,295	—	1,491	—	2,786
EX
EX rank 6	521	(47.8%)	569	(52.2%)	1090
No EX	774	(45.6%)	922	(54.4%)	1696
Total	1,295	—	1,491	—	2,786

A better analysis is to split the data and look at the different combinations of these three factors that were done. Figure 2 shows for each treatment the observed proportion of great success together with binomial confidence intervals. There was only minor differences between the treatments in the rate of success (44.8–49.2%), and not more than can be expected by chance (\(\chi^{2}\) = 2.75, df = 4, p = 0.6). The combination of refined skill, Oxford skill active and EX equipped was a little higher than the others. However, this is not more than can be expected by random variation and would if real need to be explained by the presence of some unlikely interaction effect. There is therefore no evidence that any of the 3 factors influenced the rate of success.

FIGURE 2: Proportion of great success in the five experimental conditions studied, with 95% confidence intervals

3.1.2 Batches

However, the comparisons above make the assumption that each production is an independent observation. It has been suggested that the random number generator in UWO is biased and may systematically be affected by lucky/bad streaks. In that case, it may be appropriate to consider the batch as the unit to compare (especially if any differences had been found in the overall data).

Figure 3 shows the rate of success per (complete) batch, for each treatment. Overlapping data points are visualised by the size of the points (scaled by area to the number of observations). Also shown are box-plots summarising the distributions, as well as barplots of the overall proportion of great success in the treatments. The median proportion in the batches is not much different from the overall proportion in any of the treatments: as can be expected when binomial proportions are close to 0.5, there is not much skew. Again the five treatments did not differ in the rate of success, when analysed on the batch level (one-way Anova, F_4,133 = 0.835, p = 0.5).

Proportion of great success per batch of 20 Gorz, with bar plots of the overall proportion of great success in each treatment, overlaid with box plots of the distribution of batch data.

FIGURE 3: Proportion of great success per batch of 20 Gorz, with bar plots of the overall proportion of great success in each treatment, overlaid with box plots of the distribution of batch data.

3.2 Degree of success

Examples of outcomes when crafting Gorz. a) normal success, b--d) various great success with unrefined casting, and e--f) great success with refined casting

FIGURE 4: Examples of outcomes when crafting Gorz. a) normal success, b–d) various great success with unrefined casting, and e–f) great success with refined casting

Great success Gorz came in a wide variety of configurations. Some of the variation observed is illustrated in Figure 4. Figure 5 shows the overall distribution of degree of success, for the 3 parameters measured: a) attack, b) defence and c) durability.

The three distributions of great success production are very different from each other (Figure 5). The distributions for both attack and defence are right-skewed, with higher values becoming increasingly rarer, whereas durability is bimodal with peaks at the minimum and maximum values.

For attack and defence, some great success did not result in an improved value, while durability was always increased. When attack did increase, it was always to at least 52 – the value 51 never occurred – up to the maximum of 97 attack, which was obtained twice.

Some of the curious features of these distributions are explained by treatment-effects, discussed below.

FIGURE 5: Distribution of values for a) attack, b) defence and c) durabililty. Note that the scales differ among the histograms.

The distribution of improved values for each experimental group is presented in Figure 6.

Note that the histogram bin-width for attack was set at 1 in Figure 5 above to reveal more detail, but is set to 5 in figure 6a–c below so as to not obscure the general form of the distributions. Similarly, bin-width for durability was increased from 1 to 2. The details can instead more easily be seen from the plots in Figure 6d–f.

Comparing the distribution of improved values among the treatment groups in Figure 6 immediately reveal that there is a large effect of having refined the casting skill. This is the case for both durability, attack, and defence value.

In contrast, Figure 6 also clearly shows that there is no effect of having the Oxford item improvement tech nor the EX rank 6.

Distribution of values for a) attack, b) defence and c) durabililty of great success Gorz, split by treatment group. The same data is shown in d)--f) but instead of histograms, now as box-plots overlaid by the individual data points, scaled by number of observations

FIGURE 6: Distribution of values for a) attack, b) defence and c) durabililty of great success Gorz, split by treatment group. The same data is shown in d)–f) but instead of histograms, now as box-plots overlaid by the individual data points, scaled by number of observations

The average value increased by about 10 attack, 2 defence and 8 durability from having refined the skill when obtaining a great success, compared to not having refined the skill.

The maximum value increased from 70 to 97 for attack and from 14 to 19 for defence. These maximum values were rare. Durability on the other hand saw no increase in the maximum value (70), but instead the maximum changed from being rare to being the most common value.

Not only the maximum, but also the minimum improvement values changed when refining the skill. This was the case for all three variables: the minimum attack for a great success Gorz increased from 50 to 52, defence from 10 to 11, and durability from 55 to 60.

The reason for the strange bimodal distribution of durability in Figure 5 now becomes clear: a mixture of the extremely right-skewed distribution of unrefined production and the extremely left-skewed distribution of refined production. With unrefined skill, a great success yields the minimum value about half the time, whereas with refined skill it is the maximum value instead that is obtained about half the time. The other values are more or less uniformly (evenly) distributed.

Also revealed is that the attack distribution was a mixture of the right-skewed unrefined and the more symmetrical refined. Less spectacular is the effect of refinement on the defence distribution: it is shifted to the right and with a longer tail.

FIGURE 7: Kernel density estimates for the distribution of attack, defence and durability values among the five treatment groups.

Yet another way to visualize the main differences (and similarities) between the treatment groups is offered in Figure 7, where the distributions are displayed as kernel density estimates. This is however a bit misleading since the data do not come from a continuous smooth distribution. Nevertheless, it highlights the main finding: that refining the skill has a large effect on the degree of success, and Oxford and EX do not.

Source of Variation	Attack				Defence				Durability
TABLE 3: ANOVA tables from linear models of the effect of refinement, Oxford and EX on attack, defence and durability values of great success Gorz.
Source of Variation	df	Sum Sq	F	p	df	Sum Sq	F	p	df	Sum Sq	F	p
Refined	1	39417.22	545.47	0.000	1	1134.58	565.11	0.000	1	20184.79	1128.23	0.000
Oxford	1	32.34	0.45	0.504	1	1.72	0.86	0.355	1	12.89	0.72	0.396
EX	1	161.62	2.24	0.135	1	4.46	2.22	0.136	1	2.89	0.16	0.688
Residuals	1291	93291.99	—	—	1291	2591.98	—	—	1291	23096.96	—	—

Statistical modelling (summarised in Table 3) confirmed the conclusions from the graphical analyses, that there is a large effect of refinement but no effect of Oxford or EX on the degree of great success when making Gorz.

3.3 Relation between attack, defence and durability

Scatterplots of the relationship between the measured variables of great success Gorz, a) overall and b) split in refined and unrefined skill. The inferred parameter space is indicated by rectangles.

FIGURE 8: Scatterplots of the relationship between the measured variables of great success Gorz, a) overall and b) split in refined and unrefined skill. The inferred parameter space is indicated by rectangles.

At first glance, it appears to be positive relationships between the attack, defence and durability value of great success Gorz (Figure 8a). However, this is entirely an artefact of mixing two different situations with different parameter space (Figure 8b), namely production with unrefined and with refined casting skill. When this is taken into account there is no relationship between attack, defence and durability (or combinations of these, see below).

The inferred parameter space is indicated by rectangles in the plots that are split by refined/unrefined. Note that this wasn’t always realised (for example, no item achieved both maximum attack and defence values).

Looking at each measured parameter against the sum of the other two parameters (Figure 9) revealed again that any relationship disappeared entirely when refinement was taken into account.

Figures 8 and 9 show therefore that the improvable parameters are improved independently and are not traded off against each other, and do not reflect some underlying general tier of great success.

Scatterplots of the realtionship between combinations of the measured variables and the remaining variable, a) overall and b) split in refined and unrefined skill.

FIGURE 9: Scatterplots of the realtionship between combinations of the measured variables and the remaining variable, a) overall and b) split in refined and unrefined skill.

Scatter-plots between the parameters for all five treatments, not only refined/unrefined, are shown in Figure 10 below. For ease of comparison, all the data irrespective of treatment is displayed in the background in grey. Figure 10a shows attack–defence, 10b defence–durability, 10c attack–durability.

Scatterplots of the relationship between the measured variables of great success Gorz for each treatment, a) attack--defence, b) durability--defence, c) attack--durability. Data for all treatments is displayed in the background.

FIGURE 10: Scatterplots of the relationship between the measured variables of great success Gorz for each treatment, a) attack–defence, b) durability–defence, c) attack–durability. Data for all treatments is displayed in the background.

3.4 Variation in success among batches compared to expectations

Distribution of the number of great success per batch of 20 trials. In a) is shown the observed frequency in the 138 bathces made, and b) shows one random simulation, given the observed probability of success and 138 batches of 20 trials each.

FIGURE 11: Distribution of the number of great success per batch of 20 trials. In a) is shown the observed frequency in the 138 bathces made, and b) shows one random simulation, given the observed probability of success and 138 batches of 20 trials each.

We can also ask if the variation in success rate among batches is as expected from an unbiased random process. For example, if some batches are “good” and some are “bad” to a larger extent than is likely to be observed if each trial is independently drawn from the same distribution. In that case, we should find larger variation among batches than expected from a purely random process.

The random distribution of success/failure (called a Bernoulli trial) is described by the binomial distribution, for a given probability of success and number of trials (here: batch size). The observed distribution of success in the batches is shown in Figure 11a, together with the binomial probability mass distribution.¹

The theoretical expectation for the standard deviation (SD) of a binomial distribution is \(\sqrt{(np)(1-p)}\), where n is the number of trials, p is the probability of success and 1-p is the probability of failure. If the outcomes were independently drawn at random from the same binomial distribution we can expect the observed variation among batches to be close to this theoretical expectation. The observed SD of the 128 Gorz batches was 1.92, and the expected SD is 2.23.

In order to determine how close is close enough, I performed a Monte Carlo simulation. 20 observations of success/failure were drawn randomly from a binomial distribution specified by the observed success rate in the Gorz data (0.4656), and this was repeated 128 times. Then the standard deviation of these 128 observations was calculated. An example of the outcome of one such simulation is shown in Figure 11b.

Finally, this was replicated 100,000 times to obtain a simulated distribution of SD, and a p-value was then calculated as the proportion of all simulations with a SD more extreme or equal to the observed SD, multiplied with 2 to get a two-tailed p-value. In R code:

simulsd <- replicate(10^5, sd(rbinom(138, 20, 0.4656)))
p_value <- length(which(simulsd <= obs_sd))/length(simulsd[,1])*2

The results of this simulation is shown in Figure 12. The results indicate that the observed variation among batches was somewhat less than expected by a random process (two-tailed p = 0.017), and not more as the lucky/bad streak theory predicts. The observed batches were a little too uniform to have independently been drawn at random from a binomial distribution.²

Distribution of the standard deviation from 100,000 simulations of 138 batches of 20 trials each, given the observed probability of success in the Gorz data.

FIGURE 12: Distribution of the standard deviation from 100,000 simulations of 138 batches of 20 trials each, given the observed probability of success in the Gorz data.

4 Discussion

The Gorz has a remarkably wide range of improvement possible from great success, as can be seen in Figure 6.

Before refining became available, the max improvement value was +40% of the base. Thus, for the Gorz, the maximum attainable values from great success with unrefined casting is 70 durability, 70 attack and 14 defence. With refined skill, this cap seems to be unchanged for durability but is much higher for attack and defence. If we take at face value the highest stats found in this study (97 attack and 19 defence) and assume that it is rounded down to the nearest integer, the improvement cap for refined skill appears to be at 95% for both attack and defence. However, 97 attack was very rare and it can not be excluded that even higher is possible.

To confirm that there is a general +95% cap when crafting with a refined skill for both attack and defence, similar studies should be conducted for items with higher base defence stats than the Gorz, ideally with some attack too. Preferably also items crafted with with the sewing or handicraft skills. It is likely however that there is an absolute value cap of 100 that will override the % cap.

Not only the max value cap is raised by refining, but also the minimum improvement value. Whereas unrefined great success sometimes yielded an item with no improved attack and/or defence at all (durability was always improved), refined great success always resulted in all 3 stats being improved.

The way the effect of refining is implemented is rather complex as it varies somewhat among the improvable parameters. I am a little impressed by this.

Perhaps surprisingly, there was no effect of EX or Oxford item skill on the degree or rate of great success.

This is likely because the effect of these concern only recipes that improve an already existing item, such as when forging boots, gloves or weapons. This is contrary to claims made by Papaya on behalf of Koei. To speculate, a possible resolution of the conflicting information is that the special production EX does indeed increase the rate of success but only when forging, while the Oxford skill increases the degree of improvement when forging. Thus, a practical consequence could be that the special production EX serves no purpose when using the Master’s Secrets consumable (100% success), which is routinely used when forging with expensive materials and/or items where failure results in loss of the item. More research is clearly needed.

5 Acknowledgements

6 Appendix

6.1 R code for this report

knitr::opts_chunk$set(fig.width = 5, fig.align = 'center', 
                      collapse = TRUE, warning = FALSE, message = FALSE)
#packages
library(ggplot2)
library(dplyr)
library(tidyr)
library(cowplot)
library(ggh4x)
library(janitor)
library(kableExtra)
library(fitdistrplus)
library(purrr)
library(scales)
library(gt)
library(tidyselect)
library(stringr)

#functions
#function to check for and place legend in empty facet if available
shift_legend3 <- function(p) {
  pnls <- cowplot::plot_to_gtable(p) |> gtable::gtable_filter("panel") |>
    with(setNames(grobs, layout$name)) |> purrr::keep(~identical(.x, zeroGrob()))
  
  if(length(pnls) == 0) stop("No empty facets in the plot")
  
  lemon::reposition_legend(p, "center", panel=names(pnls))
}

#ggplot settings
update_geom_defaults("bar", list(fill = I("#377eb8"), colour =("grey30")))
theme_update(plot.background = element_rect(fill = "grey92", 
                                            colour = NA), 
             #linewidth = 0.5),
             #legend.background=element_blank()
             legend.background = element_rect(fill = "grey92", colour = NA),
             legend.title = element_text(face = "bold"),
             legend.key = element_rect(color = "white", fill = "grey97")
)

# scale_y_continuous(expand = expansion(mult = c(0, .05)))

#main data
gorz <- readODS::read_ods("gorz.ods", sheet =1)
gorz$Ex_rank <- recode(as.factor(gorz$Ex_rank), "0" = "no Ex", "6" = "EX r6")
gorz$oxford_item <- recode(as.factor(gorz$oxford_item), "no" = "no Oxf", "yes" = "Oxford")
gorz$refined <- recode(as.factor(gorz$refined), "no" = "unrefined", "yes" = "refined")
gorz$great_success <- factor(gorz$great_success, levels = c("yes", "no"))
gorz$great_success2 <- factor(gorz$great_success, levels = c("no", "yes"))
gorz$condition <- with(gorz,  interaction(refined, oxford_item, Ex_rank, lex.order = TRUE, sep = ":"))
gorz <- gorz |> replace_na(list(durability = 50, attack = 50, defence =10))

treatmentTable <- data.frame(
   stringsAsFactors = FALSE,
                            treatment = c("A","B","C","D","E",
                                          "F (not done)","G (not done)",
                                          "H (not done)"),
                              casting = c("Refined","Refined","Refined",
                                          "Unrefined","Unrefined",
                                          "Unrefined","Unrefined","Refined"),
                               Oxford = c("Yes","Yes","No","Yes",
                                          "Yes","No","No","No"),
                                   EX = c("Yes","No","No",
                                          "Yes","No","Yes",
                                          "No","Yes"),
                         greatSuccess = c(261L,260L,260L,260L,254L,
                                          0L,0L,0L),
                        normalSuccess = c(269L,316L,320L,300L,286L,
                                          0L,0L,0L),
    completeBatches = c(26L, 28L, 29L, 28L, 27L, 0L, 0L, 0L),
  incompleteBatches = c(1L, 1L, 0L, 0L, 0L, 0L, 0L, 0L)
                  )
treatmentTable <- treatmentTable |>
  mutate(TotalProduced = greatSuccess+normalSuccess)

# TableTreatment2 <- kbl(treatmentTable[,-c(5,6)],
#   label = "TableTreatment",
#   caption = 'Experimental treatment groups and sample sizes',
#   col.names = c("treatment", "casting", "Oxford", "EX", "complete batches", "incomplete batches", "total produced")
# ) |> 
#   column_spec(c(5:7), width = "2em") |>
#   kable_styling(bootstrap_options = c("hover", "striped", "condensed"), full_width = F)
knitr::include_graphics('gorz recipe.png')

actual_colnames <- colnames(treatmentTable)
desired_colnames <- actual_colnames |> 
  str_remove("Batches") |> 
  str_remove("Success") |>
  str_to_title()

names(desired_colnames) <- actual_colnames

TableTreatment <- treatmentTable |>
  gt(rowname_col = "treatment") |> 
  cols_label(.list = desired_colnames) |>
  cols_label(TotalProduced = "Total Produced") |>
    tab_spanner(
    label = md('**Success**'),
    columns = 5:6) |> 
  tab_spanner(
    label = md('**Batches**'),
    columns = 7:8) |>
    tab_spanner(
    label = md('**Factors**'),
    columns = 2:4) |>
  sub_zero(zero_text = "-") |>
    grand_summary_rows(
    columns = c(5:9),
    fns = list(Total ~ sum(.))
  ) |>
    tab_options(
    data_row.padding = px(2),
    summary_row.padding = px(3), # A bit more padding for summaries
    row_group.padding = px(4)    # And even more for our groups
  ) |> 
  tab_stubhead(label = "Treatment") |>
  cols_hide(columns = c("greatSuccess", "normalSuccess")) |>
  tab_header(title = md("*TABLE 1: Experimental treatment groups and sample sizes*")) |>
  tab_options(heading.title.font.size = px(16)) |>
    tab_style(
    style = cell_text(color = "#777"),
    locations = cells_title ()) |>
  opt_stylize(style = 3, color = 'blue') |>
   tab_style(
    locations = cells_title (),
    style = cell_borders(style = 'hidden'))
#There are six color variations: "blue", "cyan", "pink", "green", "red", and "gray".
  #  tab_style(
  #   locations = cells_column_spanners(spanners = everything()),
  #   style = cell_borders(
  #       sides= c("bottom"), color = 'white', weight = px(2))
  # )
  #Bug i GT! see https://github.com/rstudio/gt/issues/648

TableTreatment
#gorz$Great_Success <- factor(gorz$great_success, levels = c("yes", "no"))

TableSuccess <- gorz |>
  rename(refinement = refined, Oxford = oxford_item) |>
  #mutate(Condition = interaction(Refinement, Oxford, Ex_rank)) |>
  tabyl(great_success, condition) |> 
  adorn_totals(c("row","col"))|>
  adorn_percentages("col") |>
  adorn_pct_formatting(rounding = "half up", digits = 1) |>
  adorn_ns(position = "front") |>
  adorn_title(placement = "top") |>
  kbl()

TableRefined <- gorz |>
  rename(casting =refined) |>
  tabyl(great_success,  casting) |> 
  adorn_totals(c("row"))|>
  adorn_percentages("col") |>
  adorn_pct_formatting(rounding = "half up", digits = 1) |>
  adorn_ns(position = "front") |>
  adorn_title(placement = "top")

TableOxford <- gorz |>
  rename(Oxford = oxford_item) |>
  tabyl(great_success, Oxford) |> 
  adorn_totals(c("row"))|>
  adorn_percentages("col") |>
  adorn_pct_formatting(rounding = "half up", digits = 1) |>
  adorn_ns(position = "front") |>
  adorn_title(placement = "top")

TableEX <- gorz |>
  rename(EX = Ex_rank) |>
  tabyl(great_success,  EX) |> 
  adorn_totals(c("row", "col"))|>
  adorn_percentages("col") |>
  adorn_pct_formatting(rounding = "half up", digits = 1) |>
  adorn_ns(position = "front") |>
  adorn_title(placement = "top")

TableEX <- TableEX[,-1]
colnames(TableEX)[3] = ""

Tabletot <- kbl(
  list(TableRefined, TableOxford[,-1], TableEX),
  label = "Tabletot",
  caption = 'Proportion great success. Note that the the totals are the same in all the 3 comparisons.',
  valign = 't'
) |>
  #column_spec(c(1), width = "8cm") |> #, width = "2em" if full_width = F to force % under the numbers
  kable_styling(bootstrap_options = c("hover", "striped"), full_width = T) 

#Redo all tables in GT once and for all

#stat tests of proportions
propTestRefined <- prop.test(x=c(514, 781), n=c(1100, 1686))
propTestOxford <- prop.test(x=c(260, 1035), n=c(580, 2206))
propTestEx <- prop.test(x=c(774, 521), n=c(1696, 1090))

gorzYesNo <- gorz |>
  group_by(condition) |>
  summarise(success = sum(great_success == "yes"),
            failure = sum(great_success == "no"))

M <- as.table(rbind(c(254, 286), c(260, 300), c(260, 320), c(260, 316), c(261, 269)))
dimnames(M) <- list(treatment = c("A", "B" ,"C", "D", "E"),
                    success = c("Yes", "No"))
chisq.test(M)

propTestRes <- prop.test(x=c(254, 260, 260, 260, 261), n=c(540, 560, 580, 576, 530))
# new Tabletot in GT

successTable <- gorz |>
  group_by(great_success) |>
  summarise(Refined = sum(refined == "refined"),
            Unrefined = sum(refined == "unrefined"),
            Oxford = sum(oxford_item == "Oxford"),
            No_oxf = sum(oxford_item == "no Oxf"),
            Ex_r6 = sum(Ex_rank == "EX r6"),
            No_EX = sum(Ex_rank == "no Ex")) |>
  pivot_longer(cols = c(-great_success), names_to = "Group") |> 
  pivot_wider(names_from = great_success) |>
  mutate(perc_yes = sprintf("(%.1f%%)", yes/(yes+no)*100)) |>
  mutate(perc_no = sprintf("(%.1f%%)", no/(yes+no)*100)) |>
  mutate(Total = yes+no) 
Comparison <- c("Casting", "Casting", "Oxford", "Oxford", "EX", "EX")
Group <- c("Refined", "Unrefined", "Oxford", "No Oxford", "EX rank 6", "No EX")
successTable <-  cbind(Group, Comparison, successTable[,-1])

TableTot2 <- successTable |> gt(groupname_col = 'Comparison',
                   rowname_col = "Group") |> 
  cols_label(yes = "Great",
             no = "Normal",
             perc_yes = "",
             perc_no = "",
             Total = "Total") |>
    tab_spanner(
    label = md('**Success**'),
    columns = 3:6) |> 
    cols_move(perc_yes, yes) |>
    cols_move(perc_no, no) |>
  summary_rows(columns = c("yes", "no", "Total"),
    fns = list("Total" ~ sum(.)),
    formatter = fmt_number,
    decimals = 0
  ) |>
    tab_options(
    data_row.padding = px(2),
    summary_row.padding = px(3), # A bit more padding for summaries
    row_group.padding = px(4),    # And even more for our groups
    row_group.font.weight = "bold"
  ) |> 
  tab_stubhead(label = "Comparison") |>
  tab_header(title = md("*TABLE 2: Proportion great success. Note that the the totals are the same in all the 3 comparisons*")) |>
  tab_options(heading.title.font.size = px(16),
              table.width = pct(60)) |>
    tab_style(
    style = cell_text(color = "#777"),
    locations = cells_title ()) |>
  opt_stylize(style = 3, color = 'blue') |>
   tab_style(
    locations = cells_title (),
    style = cell_borders(style = 'hidden'))


# ---------errorbar plot

gorzProp <- gorz |>
  group_by(condition) |>
  summarise(success = sum(great_success == "yes"),
            total = n())

#binomial confidence intervals
binomCI <- binom::binom.wilson(gorzProp$success, gorzProp$total)

CI <- cbind(gorzProp, binomCI)

propsucsCI <- CI |> 
  ggplot(aes(x = mean, y = method, colour = condition)) +
  ylab("experimental treatment") +
  xlab("proportion great success +/- 95% binomial C.I.") +
  geom_pointrange(aes(xmin = lower, xmax = upper), 
                  size =1, linewidth = 1.5) +
  scale_colour_brewer(palette = "Set1", guide = "none") +
  scale_x_continuous(breaks = seq(0, 1, 0.02), limits = c(0.4, 0.54))+ #elelr limits (0,1)
  scale_y_discrete(breaks = NULL)+
  facet_wrap(vars(condition), strip.position = "top", ncol = 1, scales = "fixed") +
  geom_text(aes(label = paste0(round(mean*100, 1), "%"), vjust = 1.5, hjust = -1/4)) +
  geom_text(aes(label = paste0(x,"/",n), vjust = -0.5, hjust = -1/4))

TableTot2
propsucsCI
batches <- gorz |>
  group_by(batch, condition) |>
  summarise(tot = n(), gs= sum(great_success == "yes")) |>
  mutate(prop = gs/tot)

pdot <- ggplot(batches, aes(condition, prop))+
  geom_count(alpha=0.9) +
  scale_size_area(n.breaks = 6)+
  guides(x = "axis_nested")

summaryBatches <- batches |>
  group_by(condition, tot) |>
  summarise(n = n())
batches2 <- batches |>
  filter(tot == 20)
summary(batches2)

gorzCompleteBatches <- gorz |>
  filter(batch != 27 & batch != 52)

baranddots <- ggplot()+
  geom_bar(data = gorzCompleteBatches, aes(x = condition,
                                           fill = factor(great_success, levels = c("no", "yes"))), 
           position = "fill", alpha = 1) +
  geom_boxplot(data = batches2,
               mapping = aes(x = condition, y = prop), outliers = F, width = 0.5, staplewidth = 0.7,  alpha =0.4)+
  geom_count(data = batches2, mapping = aes(x = condition, y=prop), alpha = 0.7) +
  scale_size_area(n.breaks = 6) +
  ylab("proportion")+
  xlab("experimental treatment")+
  scale_fill_brewer(palette = "Paired") +
  #labs(size = "n batches")+
  scale_y_continuous(limits = c(0, 1),
                     expand = expansion(mult = c(0, 0.05)))+
  guides(x = guide_axis_nested(delim = ":", inv = TRUE),
         fill = guide_legend(title = "great success", order = 1),
         size = guide_legend(title = "n batches", order = 2)) +
  theme(axis.ticks.x = element_blank()) 

# difference in success per batch
kruskal.test(batches2$prop, batches2$condition)
anova(lm(prop~condition, data=batches2))
baranddots
picA <- ggdraw() +
  draw_image("gorz1.png", scale = 0.95) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
picB <- ggdraw() +
  draw_image("gorz2.png", scale = 0.90) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
picC <- ggdraw() +
  draw_image("gorz3.png", scale = 0.95) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
picD <- ggdraw() +
  draw_image("gorz4.png", scale = 0.95) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
picE <- ggdraw() +
  draw_image("gorz19def.png", scale = 0.95) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
picF <- ggdraw() +
  draw_image("gorz5.png", scale = 0.95) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

gridpic <- plot_grid(picA, picB, picC, picD, picE, picF,#legend, 
                     ncol = 3, 
                     labels = "AUTO") +
                     #rel_widths = c(1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
#knitr::include_graphics('Examples.png')
gridpic
gshist1 <- gorz |> ggplot(aes(attack, fill = great_success2)) +
  geom_histogram(binwidth = 1) +
  scale_x_continuous(breaks = seq(50, 90, 10)) +
  facet_wrap(vars(great_success2), scales = "free_y") +
  #scale_x_continuous(sec.axis = sec_axis(~ . , name = "great success", breaks = NULL, labels = NULL))+
  scale_fill_manual(values = c("#fc9272", "#de2d26"), name = "great success") +
  theme(strip.text = element_blank(),
        legend.position = "inside",
        legend.position.inside = c(1/6, 1/3),
        legend.justification.inside = c(0, 0),
        legend.margin = NULL)

gshist2 <- gorz |> ggplot(aes(defence, fill = great_success2)) +
  geom_histogram(binwidth = 1)+
  scale_x_continuous(breaks = seq(10, 18, 2))+
  facet_wrap(vars(great_success2), , scales = "free_y")+
  #scale_x_continuous(sec.axis = sec_axis(~ . , name = "great success", breaks = NULL, labels = NULL))+
  scale_fill_manual(values = c("#9ecae1", "#3182bd"))+
  theme(strip.text = element_blank(),
        legend.position = "inside",
        legend.position.inside = c(1/6, 1/3),
        legend.justification.inside = c(0, 0),
        legend.margin = NULL) +
  guides(fill = guide_legend(title = NULL))

gshist3 <- gorz |> ggplot(aes(durability, fill = great_success2)) +
  geom_histogram(binwidth = 1)+
  scale_x_continuous(breaks = seq(50, 70, 4)) + 
  facet_wrap(vars(great_success2), scales = "free_y")+
  #scale_x_continuous(sec.axis = sec_axis(~ . , name = "great success", breaks = NULL, labels = NULL))+
  scale_fill_manual(values = c("#a1d99b", "#31a354"))+
  theme(strip.text = element_blank(),
        legend.position = "inside",
        legend.position.inside = c(1/6, 1/3),
        legend.justification.inside = c(0, 0),
        legend.margin = NULL) +
  guides(fill = guide_legend(title = NULL))

gridgshist <- plot_grid(gshist1, gshist2, gshist3,#legend, 
                   ncol = 1, 
                   labels = "AUTO", 
                   rel_heights = c(1, 1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
#for A5 portrait orentation

# ---------------------
gorz2 <- gorz[gorz$great_success=="yes",]

ItemLabels <- gorz2 |>
  group_by(refined, Ex_rank, oxford_item) |>
  summarise(n = n(), snittAtt = mean(attack), snittDef = mean(defence), snittDura = mean(durability))
ItemLabels$n2 <- paste("n =", ItemLabels$n)
ItemLabels$snittAtt2 <- sprintf("bar(x) == %.1f", ItemLabels$snittAtt)
ItemLabels$snittDef2 <- sprintf("bar(x) == %.1f", ItemLabels$snittDef)
ItemLabels$snittDura2 <- sprintf("bar(x) == %.1f", ItemLabels$snittDura)
ItemLabels$snittAndNAtt <- sprintf(
  "bar(x) == %.1f * ',' ~~ n == %.0f",
  ItemLabels$snittAtt,
  ItemLabels$n
)
#or, just make two labels!

h1 <- ggplot(gorz2, aes(attack)) +
  geom_histogram(binwidth = 5, aes(fill = condition)) +
  scale_x_continuous(breaks = seq(50, 90, 10)) + 
  facet_wrap(vars(refined, Ex_rank, oxford_item), ncol = 1, strip.position = "right") +
  scale_fill_brewer(palette = "Reds", guide = "none") +
  #theme(strip.text.y.right = element_text(angle = 0)) +
  geom_label(x = 80, y = 90, aes(label = snittAndNAtt), data = ItemLabels, parse = T, 
           hjust = "center", size = 3, label.size = NA) +
  theme(strip.text.y = element_blank())


h2 <- ggplot(gorz2, aes(defence)) +
  geom_histogram(binwidth = 1, aes(fill = condition)) +
  scale_x_continuous(breaks = seq(10, 18, 2)) + 
  facet_wrap(vars(refined, Ex_rank, oxford_item), ncol = 1, strip.position = "right") +
  scale_fill_brewer(palette = "Blues", guide = "none") +
  geom_label(x = 16, y = 100, aes(label = snittDef2), data = ItemLabels, parse = T, 
             hjust = "center", size = 3, label.size = NA) +
  theme(strip.text.y = element_blank())

h3 <- ggplot(gorz2, aes(durability)) +
  geom_histogram(binwidth = 2, aes(fill = condition)) +
  scale_x_continuous(breaks = seq(54, 70, 2)) + 
  facet_wrap(vars(refined, Ex_rank, oxford_item), ncol = 1, strip.position = "right") +
  scale_fill_brewer(palette = "Greens", guide = "none")+#, guide = "none") +
  geom_label(x = 64, y = 125, aes(label = snittDura2), data = ItemLabels, parse = T, 
             hjust = "center", size = 3, label.size = NA) +
  theme(strip.text.y = element_blank())

#for getting legend or facet to name the experimental conditions
h3b <- ggplot(gorz2, aes(durability)) +
  geom_histogram(binwidth = 2, aes(fill = condition)) +
  scale_x_continuous(breaks = seq(54, 70, 2)) + 
  facet_wrap(vars(refined, Ex_rank, oxford_item), ncol = 1, strip.position = "right") +
  scale_fill_brewer(palette = "Greens", guide = "none")+
  geom_label(x = 64, y = 125, aes(label = snittDura2), data = ItemLabels, parse = T, 
             hjust = "center", size = 3, label.size = NA)+
  theme(strip.text.y = element_text(margin = margin(1,1,1,1, "pt")))#+ 
  #theme(strip.text.y = element_text(angle = 0))
#theme(strip.text.y.right = element_text(angle = 0)) 
#how to rotate AND stack text on top of each other in the facet labels?

grid1 <- plot_grid(h1, h2, h3b,#legend, 
          ncol = 3, 
          labels = "AUTO", 
          rel_widths = c(1, 1, 1.25)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))
gorzA <- ggplot(gorz2, aes(x = condition, y = attack)) +
  geom_boxplot(aes(fill = condition), outliers = F, staplewidth = 0.7) +
  geom_count(alpha=0.4) +
  scale_size_area(n.breaks = 6) +
  scale_fill_brewer(palette = "Reds", guide = "none") +
  scale_x_discrete(labels = NULL) +
  xlab(label = NULL) +
  theme(axis.ticks.x = element_blank())
  #guides(fill = "none")

#legend <- get_legend(gorzA) 

gorzB <- ggplot(gorz2, aes(x = condition, y = defence)) +
  geom_boxplot(aes(fill = condition), outliers = F, staplewidth = 0.7) +
  scale_y_continuous(breaks = seq(10, 20, 2)) +
  geom_count(alpha=0.4) +
  scale_size_area(n.breaks = 6) +
  scale_fill_brewer(palette = "Blues", guide = "none") +
  scale_x_discrete(labels = NULL) +
  xlab(label = NULL) +
  #guides(fill = "none")+
  theme(axis.ticks.x = element_blank())

gorzC <- ggplot(gorz2, aes(x = condition, y = durability)) +
  geom_boxplot(aes(fill = condition), outliers = F, staplewidth = 0.7) +
  scale_y_continuous(breaks = seq(50, 70, 2)) +
  geom_count(alpha=0.4) +
  scale_size_area(n.breaks = 6) +
  scale_fill_brewer(palette = "Greens", guide = "none") +
  scale_x_discrete(labels = NULL) +
  xlab(label = NULL) +
  #guides(fill = "none")
  theme(axis.ticks.x = element_blank())

grid2 <- plot_grid(gorzA, gorzB, gorzC, NULL,#legend, 
          ncol = 4, 
          labels = c("D", "E", "F", ""), 
          rel_widths = c(1, 1, 1, 0)) + #0.25 for F, to align with histograms which have facet
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

compositeHistBox <- plot_grid(grid1, grid2,#legend, 
          ncol = 1, 
          #labels = "AUTO", 
          rel_widths = c(1, 1),
          rel_heights = c(1.5, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))


  
d1 <- ggplot(gorz2, aes(attack)) +
  #geom_histogram(binwidth = 5, aes(fill = interaction(Ex_rank, oxford_item))) +
  geom_density(aes(colour = condition), linewidth = 1, adjust = 2, trim = T) +
  scale_colour_brewer(palette = "Set1",guide = "none") +
  theme(strip.text.x = element_blank())
#scale_fill_manual(values = c("#B51F2E", "#377eb8"))

d2 <- ggplot(gorz2, aes(defence)) +
  geom_density(aes(colour = condition), linewidth = 1, adjust =2, trim = T) +
  scale_x_continuous(breaks = seq(10, 20, 2)) + 
  scale_colour_brewer(palette = "Set1", guide = "none") +
  theme(strip.text.x = element_blank())
#scale_fill_manual(values = c("#B51F2E", "#377eb8"))

d3 <- ggplot(gorz2, aes(durability)) +
  geom_density(aes(colour = condition), linewidth = 1, adjust = 2, trim = T) +
  scale_x_continuous(breaks = seq(50, 70, 2)) + 
  scale_colour_brewer(palette = "Set1")+ #, guide = "none") +
  theme(strip.text.x = element_blank()) +
  labs(colour = "experimental treatment")
#scale_fill_manual(values = c("#B51F2E", "#377eb8"))

#cowplot way, (without isolating legend, lazy!):
gridDensity1 <- plot_grid(d1, d2, d3,#legend, 
                   ncol = 3, 
                   labels = "AUTO", 
                   rel_widths = c(1, 1, 2)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

#legend over eller under istedet for på siden.

#patchwork way:
library(patchwork)
gridDensity2 <- d1+d2+d3+plot_layout(guides = "collect") & 
  theme(legend.position="top", legend.title.position = "top")
#may need to have label title on top or print as A4. or 4 x 10 inches
detach("package:patchwork", unload = TRUE)
gridgshist
compositeHistBox
gridDensity2
car::Anova(lm(attack~refined+oxford_item+Ex_rank, data = gorz2), type="III")
car::Anova(lm(defence~refined+oxford_item+Ex_rank, data = gorz2), type="III")
car::Anova(lm(durability~refined+oxford_item+Ex_rank, data = gorz2), type="III")
anovaAtt <- anova(lm(attack~refined+oxford_item+Ex_rank, data = gorz2))
anovaDef <- anova(lm(defence~refined+oxford_item+Ex_rank, data = gorz2))
anovaDur <- anova(lm(durability~refined+oxford_item+Ex_rank, data = gorz2))
kbl(anovaAtt)
kbl(anovaDef)
kbl(anovaDur)
# Convert ANOVA results into dataframes allows for easier name manipulation
a1<- data.frame(anovaAtt[,-3]) |> 
  tibble::rownames_to_column("Source of Variation") |>
  modify_at(c("Sum.Sq", "F.value"), ~round(.,2)) |> 
  modify_at(c("Pr..F."), ~crfsuite::txt_sprintf("%.3f", .)) 
a2 <- data.frame(anovaDef[,-3]) |> 
  tibble::rownames_to_column("Source of Variation") |>
  modify_at(c("Sum.Sq", "F.value"), ~round(.,2)) |> 
  modify_at(c("Pr..F."), ~crfsuite::txt_sprintf("%.3f", .)) 
a3 <- data.frame(anovaDur[,-3]) |>
  tibble::rownames_to_column("Source of Variation") |>
  modify_at(c("Sum.Sq", "F.value"), ~round(.,2)) |> 
  modify_at(c("Pr..F."), ~crfsuite::txt_sprintf("%.3f", .)) 

# Putting all into one dataframe/table
anova_results <- data.frame(
#   cbind(c("Refinement", "Oxford", "EX", "Residuals",
# "Refinement", "Oxford", "EX", "Residuals", 
# "Refinement", "Oxford", "EX", "Residuals"), 
                                  rbind(a1, a2, a3)) 
colnames(anova_results) <- c("", "df", "Sum Sq", "F", "p")
row.names(anova_results) <- NULL

# create HTML table using kableExtra
options(knitr.kable.NA = "")
anovaTable <- anova_results |> kbl("html",
                                    label = "anovaTable", 
caption = "Anova tables from linear modelling of the effect of refinement, Oxford and EX on attack, defence and durability values of great success Gorz.") |> 
  kable_styling(bootstrap_options = c("hover", "striped", "condensed"), full_width = F)|> 
  pack_rows("Attack", 1, 4, label_row_css = "background-color: #666; color: #fff;") |> # groups rows with label
  pack_rows("Defense", 5, 8, label_row_css = "background-color: #666; color: #fff;") |># groups rows with label
  pack_rows("Durability", 9, 12, label_row_css = "background-color: #666; color: #fff;") # groups rows with label
  
anova_results2 <- data.frame(rbind(c("Source of Variation", "df", "Sum Sq", "F", "p", 
                                     "df", "Sum Sq", "F", "p",
                                     "df", "Sum Sq", "F", "p"), 
                                   cbind(a1, a2[-1], a3[-1])))
colnames(anova_results2) <- c("", "", "", "","", "", "", "", "", "", "", "", "")
#row.names(anova_results2)[1] <- ""

anovaTable2 <- anova_results2 |> kbl("html", 
                                     label = "anovaTable2",
caption = "Anova tables from linear models of the effect of refinement, Oxford and EX on attack, defence and durability values of great success Gorz.")|> 
  column_spec(c(1), width = "2em") |>
  kable_styling(bootstrap_options = c("hover", "striped", "condensed"), full_width = F) |>
  add_header_above(c("", "Attack" = 4, "Defence" = 4, "Durability" = 4))

anova_results3 <- anova_results2
colnames(anova_results3) <- anova_results2[c(1),]
anova_results3 <- anova_results3[-c(1),]
anova_results3 <- data.frame(anova_results3)
gt(anova_results3)
Source.of.Variation <- c("Refined", "Oxford", "EX", "Residuals")
anovaTable3 <-  cbind(Source.of.Variation, anova_results3[,-1])

actual_colnames <- colnames(anovaTable3)
actual_colnames
desired_colnames <- actual_colnames |> 
  str_replace(".1|.2", "") |>
  str_replace_all("[.]", " ")
names(desired_colnames) <- actual_colnames
desired_colnames
anovaTable3 <- anovaTable3 |> gt(
                   rowname_col = "Source.of.Variation") |> 
  cols_label(.list = desired_colnames) |>   
    tab_spanner(
    label = md('**Attack**'),
    columns = 2:5) |> 
     tab_spanner(
    label = md('**Defence**'),
    columns = 6:9) |> 
     tab_spanner(
    label = md('**Durability**'),
    columns = 10:13) |> 
    tab_options(
    data_row.padding = px(2),
    summary_row.padding = px(3), # A bit more padding for summaries
    row_group.padding = px(4)    # And even more for our groups
  ) |> 
    sub_missing() |>
  tab_stubhead(label = "Source of Variation") |>
  tab_header(title = md("*TABLE 3: ANOVA tables from linear models of the effect of refinement, Oxford and EX on attack, defence and durability values of great success Gorz.*")) |>
  tab_options(heading.title.font.size = px(16)) |>
    tab_style(
    style = cell_text(color = "#777"),
    locations = cells_title ()) |>
  opt_stylize(style = 3, color = 'blue') |>
   tab_style(
    locations = cells_title (),
    style = cell_borders(style = 'hidden'))


anovaTable3
att_def <- ggplot(gorz2, aes(defence, attack)) +  
  geom_rect(xmin =10, xmax = 14, ymin = 50, ymax = 70, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  geom_rect(xmin =11, xmax = 19, ymin = 52, ymax = 97, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  #geom_smooth(aes(linetype = "Total"), method = "lm", colour = "black", linewidth= 2/3) +
  geom_smooth(method = "lm", aes(colour = refined), linewidth = 2/3) +
  geom_count(alpha = 0.6, aes(colour = refined)) +
  scale_size_area(breaks = c(10, 30, 50), guide = "none")+
  scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_colour_brewer(palette = "Set1", direction = 1, guide = "none") +
  guides(#colour = guide_legend(title = NULL, order = 2),
         linetype = guide_legend(title = NULL, order = 3))
#         size = guide_legend(position = "top"))

att_tot <- ggplot(gorz2, aes(attack, defence+durability)) +
  geom_rect(xmin = 50, xmax = 70, ymin =65, ymax = 84, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  geom_rect(xmin = 52, xmax = 97, ymin =71, ymax = 89, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  #geom_smooth(aes(linetype = "Total"), method = "lm", colour = "black", linewidth= 2/3) +
  geom_smooth(method = "lm", aes(colour = refined), linewidth = 2/3) +
  geom_count(alpha = 0.6, aes(colour = refined)) +
  scale_size_area(breaks = c(5, 15, 25), guide = "none")+
  #scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_colour_brewer(palette = "Set1", direction =1) +
  guides(colour = guide_legend(title = NULL, order = 2, position = "inside", direction = "vertical"),
         linetype = guide_legend(title = NULL, order = 3))+
 #        size = guide_legend(position = "top"))+
  theme(legend.position.inside = c(0.7, 0.15),
        legend.location = "plot",
        legend.margin = NULL)

att_dur <- ggplot(gorz2, aes(durability, attack)) +
  geom_rect(ymin = 50, ymax = 70, xmin =55, xmax = 70, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  geom_rect(ymin = 52, ymax = 97, xmin =60, xmax = 70, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  #geom_smooth(aes(linetype = "Total"), method = "lm", colour = "black", linewidth= 2/3) +
  geom_smooth(method = "lm", aes(colour = refined), linewidth = 2/3) +
  geom_count(alpha = 0.6, aes(colour = refined)) +
  scale_size_area(breaks = c(5, 20, 40), guide = "none")+
  #scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_colour_brewer(palette = "Set1", direction =1) +
  guides(colour = guide_legend(title = NULL, order = 2, position = "inside"),
         linetype = guide_legend(title = NULL, order = 3))+
#         size = guide_legend(position = "top"))+
  theme(legend.position.inside = c(0.1, 0.95),
        legend.location = "plot",
        legend.margin = NULL)

dur_tot <- ggplot(gorz2, aes(durability, attack+defence)) +
  geom_rect(xmin = 55, xmax = 70, ymin =60, ymax = 84, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  geom_rect(xmin = 60, xmax = 70, ymin =66, ymax = 116, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  #geom_smooth(aes(linetype = "Total"), method = "lm", colour = "black", linewidth= 2/3) +
  geom_smooth(method = "lm", aes(colour = refined), linewidth = 2/3) +
  geom_count(alpha = 0.6, aes(colour = refined)) +
  scale_size_area(breaks = c(5, 15, 25), guide = "none")+
  expand_limits(y = 115)+
  scale_y_continuous(breaks = seq(60, 120, 10)) +
  scale_colour_brewer(palette = "Set1", direction =1, guide = "none") 
  # guides(#colour = guide_legend(title = NULL, order = 2),
  #        linetype = guide_legend(title = NULL, order = 3),
  #        size = guide_legend(position = "top"))

def_dur<- ggplot(gorz2, aes(defence, durability)) +
  geom_rect(ymin = 55, ymax = 70, xmin =10, xmax = 14, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  geom_rect(ymin = 60, ymax = 70, xmin =11, xmax = 19, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  #geom_smooth(aes(linetype = "Total"), method = "lm", colour = "black", linewidth= 2/3) +
  geom_smooth(method = "lm", aes(colour = refined), linewidth = 2/3) +
  geom_count(alpha = 0.6, aes(colour = refined)) +
  scale_size_area(breaks = c(10, 50, 100), guide = "none")+
  scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_y_continuous(breaks = seq(50, 70, 2)) +
  scale_colour_brewer(palette = "Set1", direction =1, guide = "none") 
  # guides(#colour = guide_legend(title = NULL, order = 2, position = "top"),
  #        linetype = guide_legend(title = NULL, order = 3),
  #        size = guide_legend(position = "top"))

def_tot <- ggplot(gorz2, aes(defence, attack+durability)) +
  geom_rect(xmin = 10, xmax = 14, ymin =105, ymax = 140, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  geom_rect(xmin = 11, xmax = 19, ymin =112, ymax = 167, 
            colour = "grey30", alpha=0, linewidth = 0.5, linetype = 1)+
  #geom_smooth(aes(linetype = "Total"), method = "lm", colour = "black", linewidth= 2/3) +
  geom_smooth(method = "lm", aes(colour = refined), linewidth = 2/3) +
  geom_count(alpha = 0.6, aes(colour = refined)) +
  scale_size_area(breaks = c(5, 15, 25), guide = "none")+
  scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_y_continuous(breaks = seq(100, 170, 10)) +
  scale_colour_brewer(palette = "Set1", direction =1, guide = "none") 
  # guides(#colour = guide_legend(title = NULL, order = 2),
  #        linetype = guide_legend(title = NULL, order = 3),
  #        size = guide_legend(position = "top"))

gridRefine1 <- plot_grid(att_dur, att_def, def_dur,#legend, 
                   ncol = 3, 
                   #labels = c("B", "", ""), 
                   rel_widths = c(1, 1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

gridRefine2 <- plot_grid(def_tot, dur_tot, att_tot,#legend, 
                   ncol = 3, 
                   #labels = c("B", "", ""), 
                   rel_widths = c(1, 1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

att_def2 <- ggplot(gorz2, aes(defence, attack)) +  
  geom_smooth(method = "lm", colour = "black", linewidth= 2/3) +
  geom_count(alpha = 0.6) +
  scale_size_area(breaks = c(10, 30, 50))+
  scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_colour_brewer(palette = "Set1", direction = 1, guide = "none") +
  guides(size = guide_legend(position = "top"))

att_tot2 <- ggplot(gorz2, aes(attack, defence+durability)) +
  geom_smooth(method = "lm", colour = "black", linewidth= 2/3) +
  geom_count(alpha = 0.6) +
  scale_size_area(breaks = c(5, 15, 25))+
  guides(size = guide_legend(position = "top"))

att_dur2 <- ggplot(gorz2, aes(durability, attack)) +
  geom_smooth(method = "lm", colour = "black", linewidth= 2/3) +
  geom_count(alpha = 0.6) +
  scale_size_area(breaks = c(5, 20, 40))+
  guides(size = guide_legend(position = "top"))

dur_tot2 <- ggplot(gorz2, aes(durability, attack+defence)) +
  geom_smooth(method = "lm", colour = "black", linewidth= 2/3) +
  geom_count(alpha = 0.6) +
  scale_size_area(breaks = c(5, 15, 25))+
  expand_limits(y = 115)+
  scale_y_continuous(breaks = seq(60, 120, 10)) +
  guides(size = guide_legend(position = "top"))

def_dur2 <- ggplot(gorz2, aes(defence, durability)) +
  geom_smooth(method = "lm", colour = "black", linewidth= 2/3) +
  geom_count(alpha = 0.6) +
  scale_size_area(breaks = c(10, 50, 100))+
  scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_y_continuous(breaks = seq(50, 70, 2), limits = c(55, 70), expand = c(0.05, 0.05)) +
  guides(size = guide_legend(position = "top"))

def_tot2 <- ggplot(gorz2, aes(defence, attack+durability)) +
  geom_smooth(method = "lm", colour = "black", linewidth= 2/3) +
  geom_count(alpha = 0.6) +
  scale_size_area(breaks = c(5, 15, 25))+
  scale_x_continuous(breaks = seq(10, 20, 2)) +
  scale_y_continuous(breaks = seq(100, 170, 10)) +
  guides(size = guide_legend(position = "top"))
gridRefine3 <- plot_grid(att_dur2, att_def2, def_dur2,#legend, 
                   ncol = 3, 
                   #labels = c("A", "", ""), 
                   rel_widths = c(1, 1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

gridRefine4 <- plot_grid(def_tot2, dur_tot2, att_tot2,#legend, 
                   ncol = 3, 
                   #labels = c("A", "", ""), 
                   rel_widths = c(1, 1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

gridRefine5 <- plot_grid(gridRefine3, gridRefine1,#legend, 
                   nrow = 2, 
                   labels = c("A", "B"), 
                   rel_heights = c(1, 0.9)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

gridRefine6 <- plot_grid(gridRefine4, gridRefine2,#legend, 
                   nrow = 2, 
                   labels = c("A", "B"), 
                   rel_heights = c(1, 0.9)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

design <- matrix(c(1,2,NA,3, 4,5), nrow = 2, byrow = TRUE)
#  problemet med design og facet_manual er at empty facet ikke lenger regnes som tom fordi 
# det er NA der? shift_legend funker ihvertfall ikke
# så må ha as.table =F og recode factor for at det skal funke. Istedet plassere label manuelt grrr

scatter1b <- ggplot(gorz2, aes(attack, defence, colour = condition)) +
  geom_count(alpha = 2/3) +
  scale_size_area(guide= "none") +
  scale_y_continuous(breaks = seq(8, 20, 2),
                     expand = expansion(mult = .09))+ 
  scale_x_continuous(breaks = seq(40, 90, 10),
                     expand = expansion(mult = .09))+ 
  scale_colour_brewer(palette = "Set1", guide = "none")+
  guides(colour = guide_legend(title = "experimental treatment",
                               override.aes = list(alpha = 2/3, size = 4)))+
  facet_manual(~condition, design)+ 
  #facet_wrap(~condition, ncol = 3, as.table = FALSE) +
  gghighlight::gghighlight(unhighlighted_params = list(alpha =1/2)) +
  theme(strip.text = element_blank())


#extract the colour-only legend and save as object, or include with last plot
legendscatter <- get_legend(
  # create some space to the left of the legend t = 0, r = 0, b = 0, l = 12
  scatter1b + theme(legend.box.margin = margin(0, 0, 0, 40),
                    legend.justification = c(0, 0.7))
)

scatter1 <- ggplot(gorz2, aes(attack, defence, colour = condition)) +
  geom_count(alpha = 2/3) +
  scale_size_area() +
  scale_y_continuous(breaks = seq(8, 20, 2),
                     expand = expansion(mult = .09))+ 
  scale_x_continuous(breaks = seq(40, 90, 10),
                     expand = expansion(mult = .09))+ 
  scale_colour_brewer(palette = "Set1", guide = "none")+
  facet_manual(~condition, design)+
  gghighlight::gghighlight(unhighlighted_params = list(alpha =1/2)) +
  theme(strip.text = element_blank(),
        legend.position.inside = c(5/6, 3/4),
        legend.location = "plot")+
  guides(size = guide_legend(title.position = "top", 
                             title.hjust = 0.5, 
                             direction = "vertical",
                             position = "inside"))
  #endre retning på legend?

#scatter1 <- shift_legend3(scatter1)

scatter2 <- ggplot(gorz2, aes(attack, durability, colour = condition)) +
  geom_count(alpha = 2/3) +
  scale_size_area(breaks = c(1, 10, 20)) +
  scale_y_continuous(breaks = seq(50, 70, 5),
                     expand = expansion(mult = .15))+ 
  scale_x_continuous(breaks = seq(40, 90, 10),
                     expand = expansion(mult = .09))+ #+
  scale_colour_brewer(palette = "Set1", guide = "none")+
  facet_manual(~condition, design)+
  gghighlight::gghighlight(unhighlighted_params = list(alpha =1/2)) +
  theme(strip.text = element_blank(),
        legend.position.inside = c(5/6, 3/4),
        legend.location = "plot")+
  guides(size = guide_legend(title.position = "top", 
                             title.hjust = 0.5, 
                             direction = "vertical",
                             position = "inside"))

#scatter2 <- shift_legend3(scatter2)

scatter3 <- ggplot(gorz2, aes(durability, defence, colour = condition)) +
  geom_count(alpha = 2/3) +
  scale_size_area() +
  scale_y_continuous(breaks = seq(8, 18, 2),
                     expand = expansion(mult = .09))+ 
  scale_x_continuous(breaks = seq(50, 70, 5),
                     expand = expansion(mult = .09))+
  scale_colour_brewer(palette = "Set1", guide = "none")+
  facet_manual(~condition, design)+
  gghighlight::gghighlight(unhighlighted_params = list(alpha =1/2)) +
  theme(strip.text = element_blank(),
        legend.position.inside = c(5/6, 3/4),
        legend.location = "plot")+
  guides(size = guide_legend(title.position = "top", 
                             title.hjust=0.5, 
                             direction = "vertical",
                             position = "inside"))

#scatter3 <- shift_legend3(scatter3)

#cowplot way:
gridscatter <- plot_grid(scatter1, scatter3, scatter2, legendscatter,#legend, 
                   ncol = 2, 
                   labels = c("A", "B", "C", ""), 
                   rel_widths = c(1, 1, 1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

gridRefine5
gridRefine6
gridscatter
obs_sd <- sd(batches2$gs) #1.91727
obs_mean <- mean(batches2$gs) #9.311594
n <- length(batches2$gs) #138
size <- mean(batches2$tot) #20
prob <- mean(batches2$prop) #0.4655797

expected_sd <- sqrt((size*prob)*(1-prob)) 

set.seed(666)
df <- data.frame(x=rbinom(n, size, prob)) #x=number of successful trials
dbinomdata <- dbinom(0:20, 20, prob)
binomdata <- dbinom(0:20, 20, prob)*138
plot(dbinomdata) #probability mass
plot(binomdata) #expected counts

trials <- c(0:20)
cbind(trials, binomdata)
expbatches <- data.frame(trials, binomdata)
plot(expbatches)
expDbatches <- data.frame(trials, dbinomdata)
plot(expDbatches)

fit <- fitdist(data = batches2$gs, 
              dist = "binom", 
              method = "mle",
              fix.arg=list(size = 20))
summary(fit)
prob
plot(fit)
summary(fitdist(data = batches2$gs, 
              dist = "binom", 
              method = "mle",
              fix.arg=list(size = 18)))

binomdata18 <- dbinom(0:20, 18, 0.5173151)
plot(binomdata18)

fit2 <- fitdist(data = batches2$gs, 
              dist = "pois", 
              method = "mle")
summary(fit2)
plot(fit2)
obsBatches <- batches2 |>
  ggplot(aes(gs)) +  
  geom_bar(aes(fill= "Observed great success"), alpha = 1, width = 0.5) +
  geom_point(data = expbatches, aes(trials, y=binomdata, fill = "Binomial distribution"), shape =16, size = 2)+
  xlim(0, 20)+
  xlab("number of great success per batch of 20")+
  ylab("count")+
  annotate(geom = "text",
           label = paste0("n = ",n, " batches",
                         "\nmean = ", round(obs_mean, 2),
                         "\nobserved SD = ", round(obs_sd, 2),
                         "\nexpected SD = ", round(expected_sd, 2)),
           x = 13, 
           y = 30,
           hjust = 0,
           size = 3)+
  theme(plot.title = element_text(size = 10),
        legend.position= "inside",
        legend.position.inside = c(4/5, 0.7))+
       scale_fill_manual(name = NULL,
                        breaks = c("Observed great success", "Binomial distribution"),
                        values = c("#1F78B4", "black"))
#men hvorfor blir ikke det samme verdier med stat_thedensity????? er dette forklaringen på tilsynelatende uoverenstemmelese med binom? Nei, det er fordi size og prob ikke ble fixed i stat_theodensity, og n=16 er beste fit!!!!

oneSimul <- ggplot(data=df, aes(x)) + 
  stat_count(aes(fill = "One simulation"), width = 0.5)+
  geom_point(data=expbatches, aes(x=trials, y =binomdata, fill = "Binomial distribution"), 
             shape = 16, size =2)+
  xlab("number of great success")+
  xlim(0, 20)+
  ylim(0,33)+
  ylab("count")+
  annotate(geom = "text",
           label = paste0("example of one simulation of \nrbinom","(","number of trials = ",n ,", trial size = ",size ,", probability of success = ", round(prob, 4), ")"),
            x = 1, 
            y = 31,
            hjust = 0,
           size = 3)+
  theme(plot.title = element_text(size = 10),
        legend.position= "inside",
        legend.position.inside = c(4/5, 0.7))+
  scale_fill_manual(name = NULL,
                    breaks = c("One simulation", "Binomial distribution"),
                    values = c("#A6CEE3", "black"))

batchDistribution <- plot_grid(obsBatches, oneSimul,
                   ncol = 2, 
                   labels = "AUTO", 
                   rel_widths = c(1, 1)) +
  theme(plot.background = element_rect(fill = "grey92", colour = "grey92"))

ggplot(batches2, aes(x=gs)) +
  geom_histogram(aes(y=after_stat(density)), binwidth = 1, alpha =1, width = 0.5)+
  stat_theodensity(aes(colour = "binomial"), distri = "binom", geom = "point", size = 3, colour = "black",
                   fix.arg = list(size=20))+
  stat_theodensity(aes(colour = "binomial"), distri = "binom", geom = "point", size = 3, colour = "green",
                   fix.arg = list(size=18))+
  stat_theodensity(aes(colour = "binomial"), distri = "binom", geom = "point", size = 3, colour = "orange")+
  stat_theodensity(aes(colour = "binomial"), distri = "binom", geom = "point", size = 3, colour = "brown",
                   fix.arg = list(size=17))+
  stat_theodensity(aes(colour = "binomial"), distri = "binom", geom = "point", size = 3, colour = "blue",
                   fix.arg = list(size=21))+
  stat_theodensity(aes(colour = "normal"), distri = "norm", geom = "point", colour = "red", alpha=0.8, size = 0.5)+
  stat_theodensity(aes(colour = "poisson"), distri = "pois", geom = "point", size = 3, colour = "yellow",)+
  xlim(0, 20)+
  xlab("number of great success per batch of 20")
#plot(dbinom(0:20, 20, 0.4656))

#ggplot(expbatches, aes(x=trials, y = binomdata))+
#  geom_point(size =3)

obs_sd <- sd(batches2$gs) #1.91727
obs_mean <- mean(batches2$gs) #9.311594
size <- mean(batches2$tot) #20
n <- length(batches2$gs) #138
prob <- mean(batches2$prop) #0.4655797
simulsd <- data.frame(x = replicate(10^5, sd(rbinom(n=n, size=size, prob=prob))))
p_value <- length(which(simulsd <= obs_sd))/length(simulsd[,1])*2 #times 2 because two-tailed
simulated_sd <- mean(simulsd$x)
expected_sd <- sqrt((size*prob)*(1-prob)) 

simulSD <- ggplot(simulsd, aes(x))+
  geom_histogram(aes(fill = "simulated"), bins = 60) +
  geom_vline(aes(xintercept = obs_sd, linetype = "observed"), size =1.2, colour = "#1F78B4")+
  geom_vline(aes(xintercept = expected_sd, linetype = "theoretical"))+
  geom_vline(aes(xintercept = simulated_sd, linetype = "simulated"))+
  guides(linetype = guide_legend(title = "SD", order = 1),
         fill = guide_legend(title = NULL))+
  scale_fill_manual(values = c("#A6CEE3"))+ 
  xlab("simulated standard deviation")+
  ylab("count")+
  # labs(title = paste0("simulated distribution: replicate(10^5, sd(rbinom(138, 20, 0.4656)))", "\np (observed <= simulated) = ", round(p_value, 3)))+
  theme(plot.title = element_text(size = 10),
        legend.position= "inside",
        legend.position.inside = c(5/6, 2/3),
        legend.margin = NULL)+
  coord_cartesian(clip = 'off')+
  scale_y_continuous(labels = label_number(scale_cut = cut_short_scale()))


batchDistribution
simulsd <- replicate(10^5, sd(rbinom(138, 20, 0.4656)))
p_value <- length(which(simulsd <= obs_sd))/length(simulsd[,1])*2

simulSD

sessioninfo::session_info()
library(downloadthis)

gorz |>
  download_this(
    output_name = "GorzData",
    output_extension = ".csv",
    button_label = "Download raw data as csv",
    button_type = "primary",
    has_icon = TRUE,
    icon = "fa fa-save",
    class = "hvr-sweep-to-left"
  )

6.2 Session information

sessioninfo::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
##  setting  value
##  version  R version 4.4.1 (2024-06-14 ucrt)
##  os       Windows 10 x64 (build 19045)
##  system   x86_64, mingw32
##  ui       RTerm
##  language (EN)
##  collate  nb_NO.utf8
##  ctype    nb_NO.utf8
##  tz       Europe/Oslo
##  date     2024-11-27
##  pandoc   3.2 @ C:/Program Files/RStudio/resources/app/bin/quarto/bin/tools/ (via rmarkdown)
## 
## ─ Packages ───────────────────────────────────────────────────────────────────
##  package      * version  date (UTC) lib source
##  abind          1.4-8    2024-09-12 [1] CRAN (R 4.4.1)
##  binom          1.1-1.1  2022-05-02 [1] CRAN (R 4.3.0)
##  bookdown       0.41     2024-10-16 [1] CRAN (R 4.4.1)
##  bslib          0.8.0    2024-07-29 [1] CRAN (R 4.4.1)
##  cachem         1.1.0    2024-05-16 [1] CRAN (R 4.3.3)
##  car            3.1-3    2024-09-27 [1] CRAN (R 4.4.1)
##  carData        3.0-5    2022-01-06 [1] CRAN (R 4.3.0)
##  cellranger     1.1.0    2016-07-27 [1] CRAN (R 4.3.0)
##  cli            3.6.3    2024-06-21 [1] CRAN (R 4.4.1)
##  colorspace     2.1-1    2024-07-26 [1] CRAN (R 4.4.1)
##  commonmark     1.9.2    2024-10-04 [1] CRAN (R 4.4.1)
##  cowplot      * 1.1.3    2024-01-22 [1] CRAN (R 4.3.2)
##  crfsuite       0.4.2    2023-09-17 [1] CRAN (R 4.4.2)
##  data.table     1.16.2   2024-10-10 [1] CRAN (R 4.4.1)
##  digest         0.6.37   2024-08-19 [1] CRAN (R 4.4.1)
##  dplyr        * 1.1.4    2023-11-17 [1] CRAN (R 4.3.2)
##  evaluate       1.0.1    2024-10-10 [1] CRAN (R 4.4.1)
##  fansi          1.0.6    2023-12-08 [1] CRAN (R 4.3.2)
##  farver         2.1.2    2024-05-13 [1] CRAN (R 4.3.3)
##  fastmap        1.2.0    2024-05-15 [1] CRAN (R 4.3.3)
##  fitdistrplus * 1.2-1    2024-07-12 [1] CRAN (R 4.4.2)
##  Formula        1.2-5    2023-02-24 [1] CRAN (R 4.3.0)
##  generics       0.1.3    2022-07-05 [1] CRAN (R 4.3.0)
##  ggh4x        * 0.2.8    2024-01-23 [1] CRAN (R 4.3.2)
##  gghighlight    0.4.1    2023-12-16 [1] CRAN (R 4.3.2)
##  ggplot2      * 3.5.1    2024-04-23 [1] CRAN (R 4.3.3)
##  glue           1.8.0    2024-09-30 [1] CRAN (R 4.4.1)
##  gt           * 0.11.1   2024-10-04 [1] CRAN (R 4.4.1)
##  gtable         0.3.6    2024-10-25 [1] CRAN (R 4.4.1)
##  htmltools      0.5.8.1  2024-04-04 [1] CRAN (R 4.3.3)
##  janitor      * 2.2.0    2023-02-02 [1] CRAN (R 4.3.0)
##  jquerylib      0.1.4    2021-04-26 [1] CRAN (R 4.3.0)
##  jsonlite       1.8.9    2024-09-20 [1] CRAN (R 4.4.1)
##  kableExtra   * 1.4.0    2024-01-24 [1] CRAN (R 4.3.2)
##  knitr          1.49     2024-11-08 [1] CRAN (R 4.4.1)
##  labeling       0.4.3    2023-08-29 [1] CRAN (R 4.3.1)
##  lattice        0.22-6   2024-03-20 [1] CRAN (R 4.3.3)
##  lifecycle      1.0.4    2023-11-07 [1] CRAN (R 4.3.1)
##  lubridate      1.9.3    2023-09-27 [1] CRAN (R 4.3.1)
##  magick         2.8.5    2024-09-20 [1] CRAN (R 4.4.1)
##  magrittr       2.0.3    2022-03-30 [1] CRAN (R 4.3.0)
##  markdown       1.13     2024-06-04 [1] CRAN (R 4.3.3)
##  MASS         * 7.3-61   2024-06-13 [1] CRAN (R 4.4.1)
##  Matrix         1.7-1    2024-10-18 [1] CRAN (R 4.4.1)
##  mgcv           1.9-1    2023-12-21 [2] CRAN (R 4.4.1)
##  minty          0.0.4    2024-11-08 [1] CRAN (R 4.4.2)
##  munsell        0.5.1    2024-04-01 [1] CRAN (R 4.3.3)
##  nlme           3.1-166  2024-08-14 [1] CRAN (R 4.4.1)
##  patchwork      1.3.0    2024-09-16 [1] CRAN (R 4.4.1)
##  pillar         1.9.0    2023-03-22 [1] CRAN (R 4.3.0)
##  pkgconfig      2.0.3    2019-09-22 [1] CRAN (R 4.3.0)
##  purrr        * 1.0.2    2023-08-10 [1] CRAN (R 4.3.1)
##  R6             2.5.1    2021-08-19 [1] CRAN (R 4.3.0)
##  RColorBrewer   1.1-3    2022-04-03 [1] CRAN (R 4.3.0)
##  Rcpp           1.0.13-1 2024-11-02 [1] CRAN (R 4.4.1)
##  readODS        2.3.1    2024-11-05 [1] CRAN (R 4.4.2)
##  rlang          1.1.4    2024-06-04 [1] CRAN (R 4.3.3)
##  rmarkdown      2.29     2024-11-04 [1] CRAN (R 4.4.1)
##  rstudioapi     0.17.1   2024-10-22 [1] CRAN (R 4.4.1)
##  sass           0.4.9    2024-03-15 [1] CRAN (R 4.3.3)
##  scales       * 1.3.0    2023-11-28 [1] CRAN (R 4.3.2)
##  sessioninfo    1.2.2    2021-12-06 [1] CRAN (R 4.3.0)
##  snakecase      0.11.1   2023-08-27 [1] CRAN (R 4.3.0)
##  stringi        1.8.4    2024-05-06 [1] CRAN (R 4.3.3)
##  stringr      * 1.5.1    2023-11-14 [1] CRAN (R 4.3.1)
##  survival     * 3.7-0    2024-06-05 [1] CRAN (R 4.3.3)
##  svglite        2.1.3    2023-12-08 [1] CRAN (R 4.3.2)
##  systemfonts    1.1.0    2024-05-15 [1] CRAN (R 4.3.3)
##  tibble         3.2.1    2023-03-20 [1] CRAN (R 4.3.0)
##  tidyr        * 1.3.1    2024-01-24 [1] CRAN (R 4.3.2)
##  tidyselect   * 1.2.1    2024-03-11 [1] CRAN (R 4.3.3)
##  timechange     0.3.0    2024-01-18 [1] CRAN (R 4.3.2)
##  tzdb           0.4.0    2023-05-12 [1] CRAN (R 4.3.0)
##  utf8           1.2.4    2023-10-22 [1] CRAN (R 4.3.2)
##  vctrs          0.6.5    2023-12-01 [1] CRAN (R 4.3.2)
##  viridisLite    0.4.2    2023-05-02 [1] CRAN (R 4.3.0)
##  withr          3.0.2    2024-10-28 [1] CRAN (R 4.4.1)
##  xfun           0.49     2024-10-31 [1] CRAN (R 4.4.1)
##  xml2           1.3.6    2023-12-04 [1] CRAN (R 4.3.2)
##  yaml           2.3.10   2024-07-26 [1] CRAN (R 4.4.1)
##  zip            2.3.1    2024-01-27 [1] CRAN (R 4.3.2)
## 
##  [1] C:/R/RLibs
##  [2] C:/Program Files/R/R-4.4.1/library
## 
## ──────────────────────────────────────────────────────────────────────────────

6.3 Raw data

7 Footnotes

Converted to counts to be comparable, one could alternatively convert the histogram of observed counts to density.↩︎
However, there are a few caveats. Fitting distributions to data is both a question of choosing a probability distribution and of finding parameter estimates for that distribution. Here, I simply assumed a binomial and used the empirical proportion of success from the sample. We can use statistical optimization techniques to potentially get better estimates, but since the trial size (20) is known the maximum likelihood estimate of the probability of success will essentially be the same as the observed. Without fixing the trial size the best fit of the data was in fact a binomial distribution of 16 trials and 0.582 probability of success (AIC = 574.06) compared to AIC = 578.26 for the actual 20 trials and 0.466 probability of success, which again might suggest that something was not quite right with the UWO algorithm. We can also compare the fit with other possible discrete probability distributions. Such as hypergeometric if the trials are done without replacement, negative binomial if trials are ended when a certain number of success is achieved, geometric if trials are continued until the first success, Poisson if there is no well-defined limit of the counts, etc., but I did not pursue this further.↩︎

Crafting success in UWO – experiments with Gorz

Stinker, Anima (London), Dept. of Experimental Alchemy

Published: 23 November, 2024
Last compiled: 27 November, 2024

1 Introduction

2 Methods

3 Results

3.1 Rate of success

3.1.1 Overall and per treatment

3.1.2 Batches

3.2 Degree of success

3.3 Relation between attack, defence and durability

3.4 Variation in success among batches compared to expectations

4 Discussion

5 Acknowledgements

6 Appendix

6.1 R code for this report

6.2 Session information

6.3 Raw data

7 Footnotes

Crafting success in UWO – experiments with Gorz

Stinker, Anima (London), Dept. of Experimental Alchemy

Published: 23 November, 2024 Last compiled: 27 November, 2024

1 Introduction

2 Methods

3 Results

3.1 Rate of success

3.1.1 Overall and per treatment

3.1.2 Batches

3.2 Degree of success

3.3 Relation between attack, defence and durability

3.4 Variation in success among batches compared to expectations

4 Discussion

5 Acknowledgements

6 Appendix

6.1 R code for this report

6.2 Session information

6.3 Raw data

7 Footnotes

Published: 23 November, 2024
Last compiled: 27 November, 2024