Mercurial > hg > medcouple

annotate drawmedcouple.m @ 74:305b7361a5bd default tip @

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
showalgo: save a snapshot instead of waiting for keyboard input
author Jordi Gutiérrez Hermoso <jordigh@octave.org>
date Sun, 29 May 2016 19:05:01 -0400
parents f61da46ca3e6
children
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
36
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
1
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
2 # Set the seed, to make this "random" sampling deterministic
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
3 randg("seed", 0);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
4
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
5 z = gamrnd (3, 2, 5e3, 1);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
6
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
7 # Plot the gamma values
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
8 subplot(2,1,1);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
9 hist(z,20, "facecolor", [1, 0.5, 0.5]);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
10 title("gamma distribution with alpha=3 and beta=2", "fontsize", 20);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
11 set(gca, "fontsize", 14)
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
12
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
13 # Naive medcouple algorithm
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
14 n = length(z);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
15 n2 = ceil(n/2);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
16 z = sort(z, "descend");
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
17 zmed = median(z);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
18 z -= zmed;
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
19 zplus = z(z >= 0);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
20 zminus = z(z <= 0)';
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
21 zz = (zplus + zminus)./(zplus - zminus);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
22 zz(isnan(zz)) = 0;
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
23 zzsort = sort(zz(:), "descend");
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
24 medc_idx = ceil(n2^2/2);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
25 mc = zzsort(medc_idx);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
26
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
27 # Plot the medcouple distribution
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
28 subplot(2,1,2);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
29 hist(zzsort,20, "facecolor", [0.5, 0.5, 1]);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
30
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
31 # Draw the medcouple line on the medcouple distribution
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
32 hold on;
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
33 plot([mc,mc], [0,5e5], "linewidth", 2.5, "color", [0.95, 0.95, 0]);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
34 title("distribution of h(x_i, x_j) values", "fontsize", 20);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
35 set(gca,"xtick", [-1, -0.5, mc, 0.5, 1], "fontsize", 14);
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
36 hold off
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
37
f61da46ca3e6 drawmedcouple: new function, for Wikipedia
Jordi Gutiérrez Hermoso <jordigh@octave.org>
parents:
diff changeset
38 print medcouple-distribution.png