P-value correction: Difference between revisions

Content deleted Content added
PureFox (talk | contribs)
→‎Version 2: Changed references to Perl 6 to Raku.
Line 3,601:
Šidák
</pre>
 
=={{header|Nim}}==
{{trans|Kotlin (Version 2)}}
<lang Nim>import algorithm, math, sequtils, strformat, strutils, sugar
 
type
 
CorrectionType {.pure.} = enum
BenjaminiHochberg = "Benjamini-Hochberg"
BenjaminiYekutieli = "Benjamini-Yekutieli"
Bonferroni = "Bonferroni"
Hochberg = "Hochberg"
Holm = "Holm"
Hommel = "Hommel"
Šidák = "Šidák"
 
Direction {.pure.} = enum Up, Down
 
PValues = seq[float]
 
 
template newPValues(length: Natural): PValues =
## Create a PValues object of given length.
newSeq[float](length)
 
 
func ratchet(p: var PValues; dir: Direction) =
var m = p[0]
case dir
of Up:
for i in 1..p.high:
if p[i] > m: p[i] = m
m = p[i]
of Down:
for i in 1..p.high:
if p[i] < m: p[i] = m
m = p[i]
for i in 0..p.high:
if p[i] > 1: p[i] = 1
 
 
func schwartzian(p, mult: PValues; dir: Direction): PValues =
 
let length = p.len
let sortOrder = if dir == Up: Descending else: Ascending
let order1 = toSeq(p.pairs).sorted((x, y) => cmp(x.val, y.val), sortOrder).mapIt(it.key)
 
var pa = newPValues(length)
for i in 0..pa.high:
pa[i] = mult[i] * p[order1[i]]
 
ratchet(pa, dir)
 
let order2 = toSeq(order1.pairs).sortedByIt(it.val).mapIt(it.key)
for idx in order2:
result.add pa[idx]
 
 
proc adjust(p: PValues; ctype: CorrectionType): PValues =
let length = p.len
assert length > 0
let flength = length.toFloat
 
case ctype
 
of BenjaminiHochberg:
var mult = newPValues(length)
for i in 0..mult.high:
mult[i] = flength / (flength - i.toFloat)
return schwartzian(p, mult, Up)
 
of BenjaminiYekutieli:
var q = 0.0
for i in 1..length: q += 1 / i
var mult = newPValues(length)
for i in 0..mult.high:
mult[i] = (q * flength) / (flength - i.toFloat)
return schwartzian(p, mult, Up)
 
of Bonferroni:
result = newPValues(length)
for i in 0..result.high:
result[i] = min(p[i] * flength, 1)
return
 
of Hochberg:
var mult = newPValues(length)
for i in 0..mult.high:
mult[i] = i.toFloat + 1
return schwartzian(p, mult, Up)
 
of Holm:
var mult = newPValues(length)
for i in 0..mult.high:
mult[i] = flength - i.toFloat
return schwartzian(p, mult, Down)
 
of Hommel:
let order1 = toSeq(p.pairs).sortedByIt(it.val).mapIt(it.key)
let s = order1.mapIt(p[it])
var m = Inf
for i in 0..s.high:
m = min(m, s[i] * flength / (i + 1).toFloat)
var q, pa = repeat(m, length)
 
for j in countdown(length - 1, 2):
let lower = toSeq(0..length - j)
let upper = toSeq((length - j + 1)..<length)
var qmin = j.toFloat * s[upper[0]] / 2
for i in 1..upper.high:
let val = s[upper[i]] * j.toFloat / (i + 2).toFloat
if val < qmin: qmin = val
for idx in lower: q[idx] = min(s[idx] * j.toFloat, qmin)
for idx in upper: q[idx] = q[^j]
for i, val in q:
if pa[i] < val: pa[i] = val
 
let order2 = toSeq(order1.pairs).sortedByIt(it.val).mapIt(it.key)
return order2.mapIt(pa[it])
 
of Šidák:
result = newPValues(length)
for i in 0..result.high:
result[i] = 1 - (1 - p[i])^length
return
 
 
func pformat(p: PValues; cols = 5): string =
var lines: seq[string]
for i in countup(0, p.high, cols):
let fchunk = p[i..<(i + cols)]
var schunk = newSeq[string](fchunk.len)
for j in 0..<cols:
schunk[j] = fchunk[j].formatFloat(ffDecimal, 10)
lines.add &"[{i:2}] {schunk.join(\" \")}"
result = lines.join("\n")
 
 
func adjusted(p: PValues; ctype: CorrectionType): string =
doAssert p.len > 0 and min(p) >= 0 and max(p) <= 1, "p-values must be in range 0.0 to 1.0."
result = &"\n{ctype}\n{pformat(p.adjust(ctype))}"
 
when isMainModule:
 
const PVals = @[
4.533744e-01, 7.296024e-01, 9.936026e-02, 9.079658e-02, 1.801962e-01,
8.752257e-01, 2.922222e-01, 9.115421e-01, 4.355806e-01, 5.324867e-01,
4.926798e-01, 5.802978e-01, 3.485442e-01, 7.883130e-01, 2.729308e-01,
8.502518e-01, 4.268138e-01, 6.442008e-01, 3.030266e-01, 5.001555e-02,
3.194810e-01, 7.892933e-01, 9.991834e-01, 1.745691e-01, 9.037516e-01,
1.198578e-01, 3.966083e-01, 1.403837e-02, 7.328671e-01, 6.793476e-02,
4.040730e-03, 3.033349e-04, 1.125147e-02, 2.375072e-02, 5.818542e-04,
3.075482e-04, 8.251272e-03, 1.356534e-03, 1.360696e-02, 3.764588e-04,
1.801145e-05, 2.504456e-07, 3.310253e-02, 9.427839e-03, 8.791153e-04,
2.177831e-04, 9.693054e-04, 6.610250e-05, 2.900813e-02, 5.735490e-03]
 
for ctype in CorrectionType:
echo adjusted(PVals, ctype)</lang>
 
{{out}}
<pre style="height:60ex;overflow:scroll;">
Benjamini-Hochberg
[ 0] 0.6126681081 0.8521710465 0.1987205200 0.1891595417 0.3217789286
[ 5] 0.9301450000 0.4870370000 0.9301450000 0.6049730556 0.6826752564
[10] 0.6482628947 0.7253722500 0.5280972727 0.8769925556 0.4705703448
[15] 0.9241867391 0.6049730556 0.7856107317 0.4887525806 0.1136717045
[20] 0.4991890625 0.8769925556 0.9991834000 0.3217789286 0.9301450000
[25] 0.2304957692 0.5832475000 0.0389954722 0.8521710465 0.1476842609
[30] 0.0168363750 0.0025629017 0.0351608437 0.0625018947 0.0036365888
[35] 0.0025629017 0.0294688286 0.0061660636 0.0389954722 0.0026889914
[40] 0.0004502862 0.0000125223 0.0788155476 0.0314261300 0.0048465270
[45] 0.0025629017 0.0048465270 0.0011017083 0.0725203250 0.0220595769
 
Benjamini-Yekutieli
[ 0] 1.0000000000 1.0000000000 0.8940844244 0.8510676197 1.0000000000
[ 5] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[10] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[15] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 0.5114323399
[20] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[25] 1.0000000000 1.0000000000 0.1754486368 1.0000000000 0.6644618149
[30] 0.0757503083 0.0115310209 0.1581958559 0.2812088585 0.0163617595
[35] 0.0115310209 0.1325863108 0.0277423864 0.1754486368 0.0120983246
[40] 0.0020259303 0.0000563403 0.3546073326 0.1413926119 0.0218055202
[45] 0.0115310209 0.0218055202 0.0049568120 0.3262838334 0.0992505663
 
Bonferroni
[ 0] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[ 5] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[10] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[15] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[20] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[25] 1.0000000000 1.0000000000 0.7019185000 1.0000000000 1.0000000000
[30] 0.2020365000 0.0151667450 0.5625735000 1.0000000000 0.0290927100
[35] 0.0153774100 0.4125636000 0.0678267000 0.6803480000 0.0188229400
[40] 0.0009005725 0.0000125223 1.0000000000 0.4713919500 0.0439557650
[45] 0.0108891550 0.0484652700 0.0033051250 1.0000000000 0.2867745000
 
Hochberg
[ 0] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[ 5] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[10] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[15] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[20] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[25] 0.9991834000 0.9991834000 0.4632662100 0.9991834000 0.9991834000
[30] 0.1575884700 0.0138396690 0.3938014500 0.7600230400 0.0250197306
[35] 0.0138396690 0.3052970640 0.0542613600 0.4626366400 0.0165641872
[40] 0.0008825610 0.0000125223 0.9930759000 0.3394022040 0.0369228426
[45] 0.0102358057 0.0397415214 0.0031729200 0.8992520300 0.2179486200
 
Holm
[ 0] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[ 5] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[10] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[15] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[20] 1.0000000000 1.0000000000 1.0000000000 1.0000000000 1.0000000000
[25] 1.0000000000 1.0000000000 0.4632662100 1.0000000000 1.0000000000
[30] 0.1575884700 0.0139534054 0.3938014500 0.7600230400 0.0250197306
[35] 0.0139534054 0.3052970640 0.0542613600 0.4626366400 0.0165641872
[40] 0.0008825610 0.0000125223 0.9930759000 0.3394022040 0.0369228426
[45] 0.0102358057 0.0397415214 0.0031729200 0.8992520300 0.2179486200
 
Hommel
[ 0] 0.9991834000 0.9991834000 0.9991834000 0.9987623800 0.9991834000
[ 5] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[10] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[15] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9595180000
[20] 0.9991834000 0.9991834000 0.9991834000 0.9991834000 0.9991834000
[25] 0.9991834000 0.9991834000 0.4351894700 0.9991834000 0.9766522500
[30] 0.1414255500 0.0130434007 0.3530936533 0.6887708800 0.0238560222
[35] 0.0132245726 0.2722919760 0.0542613600 0.4218157600 0.0158112696
[40] 0.0008825610 0.0000125223 0.8743649143 0.3016908480 0.0351646120
[45] 0.0095824564 0.0387722160 0.0031729200 0.8122276400 0.1950066600
 
Šidák
[ 0] 1.0000000000 1.0000000000 0.9946598274 0.9914285749 0.9999515274
[ 5] 1.0000000000 0.9999999688 1.0000000000 1.0000000000 1.0000000000
[10] 1.0000000000 1.0000000000 0.9999999995 1.0000000000 0.9999998801
[15] 1.0000000000 1.0000000000 1.0000000000 0.9999999855 0.9231179729
[20] 0.9999999956 1.0000000000 1.0000000000 0.9999317605 1.0000000000
[25] 0.9983109511 1.0000000000 0.5068253940 1.0000000000 0.9703301333
[30] 0.1832692440 0.0150545753 0.4320729669 0.6993672225 0.0286818157
[35] 0.0152621104 0.3391808707 0.0656206307 0.4959194266 0.0186503726
[40] 0.0009001752 0.0000125222 0.8142104886 0.3772612062 0.0430222116
[45] 0.0108312558 0.0473319661 0.0032997780 0.7705015898 0.2499384839</pre>
 
=={{header|Perl}}==