## The MCPMod package
library("MCPMod")

## Preliminaries
mods1 <- list(linear = NULL, emax = 0.2, logistic = c(0.25, 0.09))
mods2 <- list(linear = NULL, emax = c(0.05, 0.2), betaMod = c(0.5, 1),
  logistic = matrix(c(0.25, 0.7, 0.09, 0.06), byrow = FALSE,
  nrow = 2))

## Planning code
guesst(d = 0.2, p = 0.9, model = "emax")
guesst(d = c(0.05, 0.2), p = c(0.2, 0.9), model = "logistic")

doses <- c(0, 0.05, 0.2, 0.6, 1)
plotModels(mods2, doses, base = 0, maxEff = 0.4, scal = 1.2)

doses <- c(0, 0.05, 0.2, 0.6, 1)
fmods2 <- fullMod(mods2, doses, base = 0, maxEff = 0.4, scal = 1.2)

pM <- planMM(mods2, doses, n = 20, alpha = 0.05, twoSide = FALSE,
  scal = 1.2)
pM
plot(pM)

pM <- powerMM(fmods2, sigma = 1, alpha = 0.05, lower = 10,
  upper = 110, step = 10)
plot(pM, line.at = 0.9, models = "none")

sampSize(fmods2, sigma = 1, sumFct = mean, power = 0.9, alpha = 0.05,
  twoSide = FALSE, upperN = 100)

mods3 <- list(linear = NULL, emax = 0.15)
Lfit <- LP(mods3, model = "emax", type = "LP1",
  paramRange = c(0.03, 0.8), len = 30, doses = doses, n = 92,
  base = 0, maxEff = 0.4, sigma = 1, alpha = 0.05, twoSide = FALSE)
plot(Lfit, spldf = 25)

## Analysis code
data("biom")
dfe <- MCPMod(biom, mods2, alpha = 0.05, dePar = 0.05, pVal = TRUE,
  selModel = "maxT", doseEst = "MED2", clinRel = 0.4, scal = 1.2)
dfe
summary(dfe)
plot(dfe, complData = TRUE, clinRel = TRUE)

dfe2 <- MCPMod(biom, mods2, alpha = 0.05, dePar = 0.95,
  selModel = "aveAIC", doseEst = "ED", scal = 1.2)
dfe2
summary(dfe2)
plot(dfe2, doseEst = TRUE)