UN Votes analysis using d3 and Numeric Javascript, forked from Principal Component Analysis by ktaneishi.
Preview on http://bl.ocks.org/Fil/2481e1f3b1d92002317e
<!DOCTYPE html>
<meta charset="utf-8">
<style>
body {
font: 11px "Cisalpin LT Std", "Lucida Grande";
}
.axis path,
.axis line {
fill: none;
stroke: #000;
shape-rendering: crispEdges;
}
.axis {display: none;}
.dot circle {
stroke: #000;
}
.nuclear {
stroke: #d24232;
stroke-width: 2px;
}
</style>
<body>
<svg>
<defs>
<filter id="blur" x="0" y="0">
<feGaussianBlur in="SourceGraphic" stdDeviation="2" />
</filter>
</defs>
</svg>
<script src="//d3js.org/d3.v3.min.js" charset="utf-8"></script>
<script src="pca.js"></script>
<script>
function parseVote(x) {
if (x == 'Y') return 0;
if (x == 'N') return 5;
if (x == 'A') return 2;
if (x == 'X') return 1;
}
var reposition = true,
zero = {
x: -2,
y: 6
},
inflation = 1.1;
var randomseed = 0,
version = 50;
Math.random = function () {
randomseed++;
return Math.sin(version * randomseed);
}
var margin = {
top: 20,
right: 20,
bottom: 30,
left: 40
};
margin = {
top: 2,
right: 2,
bottom: 2,
left: 2
};
var width = 900 - margin.left - margin.right,
height = 650 - margin.top - margin.bottom;
var x = d3.scale.linear()
.range([0, width]);
var y = d3.scale.linear()
.range([height, 0]);
var color = d3.scale.ordinal()
.domain(['AG','APG','EEG','WEOG','GRULAC'])
.range(['#c27d66', '#e3a26a','#f5ce96', '#fae7b9', '#c4b38e', 'red']);
var nuclear = function(iso){
return ['USA','RUS','FRA','GBR','PRK','CHN','PAK','IND','ISR']
.indexOf(iso)>-1;
}
var xAxis = d3.svg.axis()
.scale(x)
.orient("bottom");
var yAxis = d3.svg.axis()
.scale(y)
.orient("left");
var svg = d3.select("svg")
.attr("width", width + margin.left + margin.right)
.attr("height", height + margin.top + margin.bottom)
.append("g")
.attr("transform", "translate(" + margin.left + "," + margin.top + ")");
var groups = d3.map();
var names = d3.map();
d3.csv("country-codes.csv", function (error, data) {
data.forEach(function (d) {
groups.set(d.iso, d.group);
names.set(d.iso, d.country);
});
d3.csv("votes-session-70-nations-unies.csv", function (error, data) {
matrix = [];
data.map(function (d) {
d = d3.values(d);
d = d.slice(1).map(parseVote);
matrix.push(d);
});
var pca = new PCA();
matrix = pca.scale(matrix, true, true);
pc = pca.pca(matrix, 2);
data.map(function (d, i) {
d.x = d.pc1 = -pc[i][0];
d.y = d.pc2 = -pc[i][1];
});
svg.append('g').attr('id', 'bg');
var dots = svg.selectAll(".dot")
.data(data)
.enter().append("g")
.attr({class: "dot"
,id: function(d){return d.iso;}});
var taken = {},
step;
if (reposition) {
for (step = 0.001; step < 1; step *= inflation) {
taken = {};
dots.each(function (d) {
var angle = Math.atan((d.x - zero.x) / (d.y - zero.y));
do {
var x1 = Math.ceil(d.x / step),
y1 = Math.ceil(d.y / step),
v = "" + x1 + "-" + y1;
} while (!!taken[v] && (d.x -= (0.3 + Math.random()) * Math.sin(angle)) && (d.y -= (0.3 + Math.random()) * Math.cos(angle)))
taken[v] = 1;
});
}
step = step / inflation;
dots.each(function (d) {
d.x = step * Math.ceil(d.x / step);
d.y = step * Math.ceil(d.y / step);
});
}
dots
.attr('title',function (d) {
return names.get(d.iso);
});
dots
.append("rect")
.attr({
y: -9,
x: -9,
width: 18,
height: 12,
class: function(d){
if (nuclear(d.iso)) return 'nuclear';
}
})
.style({fill: function (d) {
return color(groups.get(d.iso));
}})
;
dots
.append("text")
.text(function (d) {
return d.iso;
})
.style({
'text-anchor': 'middle',
fill: 'black',
'font-size': '7px',
})
x.domain(d3.extent(data, function (d) {
return d.x;
}))
.nice()
.domain([-20, 17]);
y.domain(d3.extent(data, function (d) {
return d.y;
})).nice()
.domain([-22, 25]);
dots
.attr("transform", function (d) {
return "translate(" + x(d.pc1) + "," + y(d.pc2) + ")";
});
var area = d3.svg.area()
.x(function (d) {
return x(d.x);
})
.y0(function (d) {
return y(d.y);
})
.y1(0);
d3.select('#bg')
.append('rect')
.attr({
x: 0,
y: 0,
width: width,
height: height,
})
.style({
fill: '#3f151d',
stroke: '#031622',
'stroke-width': '3px'
});
d3.select('#bg')
.append("path")
.datum([{
x: -20,
y: -7
}, {
x: 12,
y: 25
}])
.attr("class", "area")
.attr("d", area)
.style({
fill: "#fef6ea"
})
.attr('filter', 'url(#blur)')
;
svg.append("g")
.attr("class", "x axis")
.attr("transform", "translate(0," + height + ")")
.call(xAxis)
.append("text")
.attr("class", "label")
.attr("x", width)
.attr("y", -6)
.style("text-anchor", "end")
.text("PC1");
svg.append("g")
.attr("class", "y axis")
.call(yAxis)
.append("text")
.attr("class", "label")
.attr("transform", "rotate(-90)")
.attr("y", 6)
.attr("dy", ".71em")
.style("text-anchor", "end")
.text("PC2")
dots
.attr("transform", function (d) {
return "translate(" + x(d.pc1) + "," + y(d.pc2) + ")";
})
//.transition()
.attr("transform", function (d) {
return "translate(" + x(d.x) + "," + y(d.y) + ")";
});
var grouplegend = svg
.append('g')
.attr('id','groups')
.attr("transform", function (d, i) {
return "translate("+(width-25)+","+ (35) + ")";
})
;
grouplegend.append('text')
.text('Groupes régionaux')
.attr({'text-anchor': 'end', fill: '#fdfef1'})
.style({'font-weight': 'bold','font-size': '13px'})
;
var legend = grouplegend
.selectAll(".legend")
.data(color.domain())
.enter().append("g")
.attr("class", "legend")
.attr("transform", function (d, i) {
return "translate(0," + (10 + i * 20) + ")";
});
var grups = {'APG': 'Asie et Pacifique', 'EEG': 'Europe de l’est', 'AG':'Afrique', 'WEOG':'Europe occidentale et autres', 'GRULAC':'Amérique latine et Caraïbes'}
legend.append("rect")
.attr("x", -20)
.attr("width", 18)
.attr("height", 18)
.style("fill", color);
legend.append("text")
.attr("x", -30)
.attr("y", 9)
.attr("dy", ".35em")
.style({
"text-anchor": "end",
fill: '#fdfef1',
})
.text(function (d) {
return grups[d] || d;
});
grouplegend.append("rect")
.attr({x: -20, y: 120, width: 18, height: 18, class:'nuclear'})
.style("fill", "#dadfbc");
grouplegend.append("text")
.attr("x", -30)
.attr("y", 129)
.attr("dy", ".35em")
.style({
"text-anchor": "end",
fill: '#fdfef1',
})
.text('Pays disposant de l’arme nucléaire');
svg.append('text')
.text('Qui vote avec qui ?')
.attr({
x: 15,
y: 40,
'font-size': '32px'
});
svg.append('g')
.selectAll('text')
.data(['Les 193 États-membres des Nations unies, répartis selon la similarité de leurs votes sur les questions de désarmement', 'et de sécurité internationale lors de la 70e Assemblée générale, septembre-décembre 2015.',
'',
'Principaux points de clivage entre les nations :',
'1. — En bas à droite : consensus, pays votant “oui” à la quasi-totalité des résolutions',
'2. — Vers la gauche : pays occidentaux, plutôt opposés aux résolutions sur le désarmement nucléaire', // L15 L38 L40 L44 L51
'3. — Vers le haut : pays plutôt opposés aux résolutions sur les armes chimiques et les bombes à sous-munitions', // L27 L26 L49
'4. — En haut à gauche : dissensus, pays votant “non” ou s’abstenant plus que les autres.',
])
.enter()
.append('text')
.attr({
transform: function(d,i) {return 'translate(15,' + (height - 118 + 14*i) + ')';
},
fill: '#fdfef1',
})
.text(function(d){return d;});
svg.append('g')
.attr('id', 'points')
.selectAll('text')
.data([
{t: '➊', x: width * 0.64, y: height * 0.7},
{t: '➋', x: width * 0.1, y: height * 0.7},
{t: '➌', x: width * 0.64, y: height * 0.13},
{t: '➍', x: width * 0.1, y: height * 0.13},
])
.enter()
.append('text')
.text(function(d) {
return d.t;
})
.attr('transform', function(d) { return 'translate(' + d.x + ',' + d.y + ')'; })
.style({
'font-size': '23px',
fill: function(d, i) {
if (i < 2) return '#fdfef1'; return '#3f151d';
},
});
config = {width: width, height:height}
graphsignature(['Philippe Rivière — Visionscarto','Décembre 2015']);
});
});
function graphsignature(lignes) {
signature = svg
.append('g');
signature
.selectAll('text')
.data(lignes)
.enter()
.append('text')
.attr("y", function (d, i) {
return 14 * i;
})
.attr("fill", "white")
.attr("text-anchor", "middle")
.style("font-size", "8px")
.text(function (d) {
return d.toUpperCase();
});
var bbox = signature[0][0].getBBox();
signature
.attr({
transform: 'translate(' + [config.width - bbox.width / 2 - 19, config.height - bbox.height -10] + ')',
});
signature
.selectAll('line')
.data(lignes.slice(1))
.enter()
.append('line')
.attr({
x1: -bbox.width / 2,
x2: bbox.width / 2,
y1: function (d, i) {
return 14 * i + 4;
},
y2: function (d, i) {
return 14 * i + 4;
}
})
.attr("stroke", "white");
}
</script>
num,country,iso,group
1,AFGHANISTAN,AFG,APG
2,ALBANIA,ALB,EEG
3,ALGERIA,DZA,AG
4,ANDORRA,AND,WEOG
5,ANGOLA,AGO,AG
6,ANTIGUA AND BARBUDA,ATG,GRULAC
7,ARGENTINA,ARG,GRULAC
8,ARMENIA,ARM,EEG
9,AUSTRALIA,AUS,WEOG
10,AUSTRIA,AUT,WEOG
11,AZERBAIJAN,AZE,EEG
12,BAHAMAS,BHS,GRULAC
13,BAHRAIN,BHR,APG
14,BANGLADESH,BGD,APG
15,BARBADOS,BRB,GRULAC
16,BELARUS,BLR,EEG
17,BELGIUM,BEL,WEOG
18,BELIZE,BLZ,GRULAC
19,BENIN,BEN,AG
20,BHUTAN,BTN,APG
21,BOLIVIA (PLURINATIONAL STATE OF),BOL,GRULAC
22,BOSNIA AND HERZEGOVINA,BIH,EEG
23,BOTSWANA,BWA,AG
24,BRAZIL,BRA,GRULAC
25,BRUNEI DARUSSALAM,BRN,APG
26,BULGARIA,BGR,EEG
27,BURKINA FASO,BFA,AG
28,BURUNDI,BDI,AG
29,CABO VERDE,CPV,AG
30,CAMBODIA,KHM,APG
31,CAMEROON,CMR,AG
32,CANADA,CAN,WEOG
33,CENTRAL AFRICAN REPUBLIC,CAF,AG
34,CHAD,TCD,AG
35,CHILE,CHL,GRULAC
36,CHINA,CHN,APG
37,COLOMBIA,COL,GRULAC
38,COMOROS,COM,AG
39,CONGO,COG,AG
40,COSTARICA,CRI,GRULAC
41,COTE D'IVOIRE,CIV,AG
42,CROATIA,HRV,EEG
43,CUBA,CUB,GRULAC
44,CYPRUS,CYP,APG
45,CZECH REPUBLIC,CZE,EEG
46,DEMOCRATIC PEOPLE'S REPUBLIC OF KOREA,PRK,APG
47,DEMOCRATIC REPUBLIC OF THE CONGO,COD,AG
48,DENMARK,DNK,WEOG
49,DJIBOUTI,DJI,AG
50,DOMINICA,DMA,GRULAC
51,DOMINICAN REPUBLIC,DOM,GRULAC
52,ECUADOR,ECU,GRULAC
53,EGYPT,EGY,AG
54,EL SALVADOR,SLV,GRULAC
55,EQUATORIAL GUINEA,GNQ,AG
56,ERITREA,ERI,AG
57,ESTONIA,EST,EEG
58,ETHIOPIA,ETH,AG
59,FIJI,FJI,APG
60,FINLAND,FIN,WEOG
61,FRANCE,FRA,WEOG
62,GABON,GAB,AG
63,GAMBIA,GMB,AG
64,GEORGIA,GEO,EEG
65,GERMANY,DEU,WEOG
66,GHANA,GHA,AG
67,GREECE,GRC,WEOG
68,GRENADA,GRD,GRULAC
69,GUATEMALA,GTM,GRULAC
70,GUINEA,GIN,AG
71,GUINEA-BISSAU,GNB,AG
72,GUYANA,GUY,GRULAC
73,HAITI,HTI,GRULAC
74,HONDURAS,HND,GRULAC
75,HUNGARY,HUN,EEG
76,ICELAND,ISL,WEOG
77,INDIA,IND,APG
78,INDONESIA,IDN,APG
79,IRAN (ISLAMIC REPUBLIC OF),IRN,APG
80,IRAQ,IRQ,APG
81,IRELAND,IRL,WEOG
82,ISRAEL,ISR,WEOG
83,ITALY,ITA,WEOG
84,JAMAICA,JAM,GRULAC
85,JAPAN,JPN,APG
86,JORDAN,JOR,APG
87,KAZAKHSTAN,KAZ,APG
88,KENYA,KEN,AG
89,KIRIBATI,KIR,APG
90,KUWAIT,KWT,APG
91,KYRGYZSTAN,KGZ,APG
92,LAO PEOPLE'S DEMOCRATIC REPUBLIC,LAO,APG
93,LATVIA,LVA,EEG
94,LEBANON,LBN,APG
95,LESOTHO,LSO,AG
96,LIBERIA,LBR,AG
97,LIBYA,LBY,AG
98,LIECHTENSTEIN,LIE,WEOG
99,LITHUANIA,LTU,EEG
100,LUXEMBOURG,LUX,WEOG
101,MADAGASCAR,MDG,AG
102,MALAWI,MWI,AG
103,MALAYSIA,MYS,APG
104,MALDIVES,MDV,APG
105,MALI,MLI,AG
106,MALTA,MLT,WEOG
107,MARSHALL ISLANDS,MHL,APG
108,MAURITANIA,MRT,AG
109,MAURITIUS,MUS,AG
110,MEXICO,MEX,GRULAC
111,MICRONESIA (FEDERATED STATES OF),FSM,APG
112,MONACO,MCO,WEOG
113,MONGOLIA,MNG,APG
114,MONTENEGRO,MNE,EEG
115,MOROCCO,MAR,AG
116,MOZAMBIQUE,MOZ,AG
117,MYANMAR,MMR,APG
118,NAMIBIA,NAM,AG
119,NAURU,NRU,APG
120,NEPAL,NPL,APG
121,NETHERLANDS,NLD,WEOG
122,NEW ZEALAND,NZL,WEOG
123,NICARAGUA,NIC,GRULAC
124,NIGER,NER,AG
125,NIGERIA,NGA,AG
126,NORWAY,NOR,WEOG
127,OMAN,OMN,APG
128,PAKISTAN,PAK,APG
129,PALAU,PLW,APG
130,PANAMA,PAN,GRULAC
131,PAPUA NEW GUINEA,PNG,APG
132,PARAGUAY,PRY,GRULAC
133,PERU,PER,GRULAC
134,PHILIPPINES,PHL,APG
135,POLAND,POL,EEG
136,PORTUGAL,PRT,WEOG
137,QATAR,QAT,APG
138,REPUBLIC OF KOREA,KOR,APG
139,REPUBLIC OF MOLDOVA,MDA,EEG
140,ROMANIA,ROU,EEG
141,RUSSIAN FEDERATION,RUS,EEG
142,RWANDA,RWA,AG
143,SAINT KITTS AND NEVIS,KNA,GRULAC
144,SAINTLUCIA,LCA,GRULAC
145,SAINT VINCENT AND THE GRENADINES,VCT,GRULAC
146,SAMOA,WSM,APG
147,SAN MARINO,SMR,WEOG
148,SAO TOME AND PRINCIPE,STP,AG
149,SAUDI ARABIA,SAU,APG
150,SENEGAL,SEN,AG
151,SERBIA,SRB,EEG
152,SEYCHELLES,SYC,AG
153,SIERRA LEONE,SLE,AG
154,SINGAPORE,SGP,APG
155,SLOVAKIA,SVK,EEG
156,SLOVENIA,SVN,EEG
157,SOLOMON ISLANDS,SLB,APG
158,SOMALIA,SOM,AG
159,SOUTH AFRICA,ZAF,AG
160,SOUTH SUDAN,SSD,AG
161,SPAIN,ESP,WEOG
162,SRI LANKA,LKA,APG
163,SUDAN,SDN,AG
164,SURINAME,SUR,GRULAC
165,SWAZILAND,SWZ,AG
166,SWEDEN,SWE,WEOG
167,SWITZERLAND,CHE,WEOG
168,SYRIAN ARAB REPUBLIC,SYR,APG
169,TAJIKISTAN,TJK,APG
170,THAILAND,THA,APG
171,THE FORMER YUGOSLAV REPUBLIC OF MACEDONIA,MKD,EEG
172,TIMOR-LESTE,TLS,APG
173,TOGO,TGO,AG
174,TONGA,TON,APG
175,TRINIDAD AND TOBAGO,TTO,GRULAC
176,TUNISIA,TUN,AG
177,TURKEY,TUR,WEOG
178,TURKMENISTAN,TKM,APG
179,TUVALU,TUV,APG
180,UGANDA,UGA,AG
181,UKRAINE,UKR,EEG
182,UNITED ARAB EMIRATES,ARE,APG
183,UNITED KINGDOM,GBR,WEOG
184,UNITED REPUBLIC OF TANZANIA,TZA,AG
185,UNITED STATES,USA,WEOG
186,URUGUAY,URY,GRULAC
187,UZBEKISTAN,UZB,APG
188,VANUATU,VUT,APG
189,VENEZUELA,VEN,GRULAC
190,VIET NAM,VNM,APG
191,YEMEN,YEM,APG
192,ZAMBIA,ZMB,AG
193,ZIMBABWE,ZWE,AG
var PCA = function(){
this.scale = scale;
this.pca = pca;
function mean(X){
// mean by col
var T = transpose(X);
return T.map(function(row){ return d3.sum(row) / X.length; });
}
function transpose(X){
return d3.range(X[0].length).map(function(i){
return X.map(function(row){ return row[i]; });
});
}
function dot(X,Y){
return X.map(function(row){
return transpose(Y).map(function(col){
return d3.sum(d3.zip(row,col).map(function(v){
return v[0]*v[1];
}));
});
});
}
function diag(X){
return d3.range(X.length).map(function(i){
return d3.range(X.length).map(function(j){ return (i == j) ? X[i] : 0; });
});
}
function zeros(i,j){
return d3.range(i).map(function(row){
return d3.range(j).map(function(){ return 0; });
});
}
function trunc(X,d){
return X.map(function(row){
return row.map(function(x){ return (x < d) ? 0 : x; });
});
}
function same(X,Y){
return d3.zip(X,Y).map(function(v){
return d3.zip(v[0],v[1]).map(function(w){ return w[0] == w[1]; });
}).map(function(row){
return row.reduce(function(x,y){ return x*y; });
}).reduce(function(x,y){ return x*y; });
}
function std(X){
var m = mean(X);
return sqrt(mean(mul(X,X)), mul(m,m));
}
function sqrt(V){
return V.map(function(x){ return Math.sqrt(x); });
}
function mul(X,Y){
return d3.zip(X,Y).map(function(v){
if (typeof(v[0]) == 'number') return v[0]*v[1];
return d3.zip(v[0],v[1]).map(function(w){ return w[0]*w[1]; });
});
}
function sub(x,y){
console.assert(x.length == y.length, 'dim(x) == dim(y)');
return d3.zip(x,y).map(function(v){
if (typeof(v[0]) == 'number') return v[0]-v[1];
else return d3.zip(v[0],v[1]).map(function(w){ return w[0]-w[1]; });
});
}
function div(x,y){
console.assert(x.length == y.length, 'dim(x) == dim(y)');
return d3.zip(x,y).map(function(v){ return v[0]/v[1]; });
}
function scale(X, center, scale){
// compatible with R scale()
if (center){
var m = mean(X);
X = X.map(function(row){ return sub(row, m); });
}
if (scale){
var s = std(X);
X = X.map(function(row){ return div(row, s); });
}
return X;
}
// translated from http://stitchpanorama.sourceforge.net/Python/svd.py
function svd(A){
var temp;
// Compute the thin SVD from G. H. Golub and C. Reinsch, Numer. Math. 14, 403-420 (1970)
var prec = Math.pow(2,-52) // assumes double prec
var tolerance = 1.e-64/prec;
var itmax = 50;
var c = 0;
var i = 0;
var j = 0;
var k = 0;
var l = 0;
var u = A.map(function(row){ return row.slice(0); });
var m = u.length;
var n = u[0].length;
console.assert(m >= n, 'Need more rows than columns');
var e = d3.range(n).map(function(){ return 0; });
var q = d3.range(n).map(function(){ return 0; });
var v = zeros(n,n);
function pythag(a,b){
a = Math.abs(a)
b = Math.abs(b)
if (a > b)
return a*Math.sqrt(1.0+(b*b/a/a))
else if (b == 0)
return a
return b*Math.sqrt(1.0+(a*a/b/b))
}
// Householder's reduction to bidiagonal form
var f = 0;
var g = 0;
var h = 0;
var x = 0;
var y = 0;
var z = 0;
var s = 0;
for (i=0; i < n; i++)
{
e[i]= g;
s= 0.0;
l= i+1;
for (j=i; j < m; j++)
s += (u[j][i]*u[j][i]);
if (s <= tolerance)
g= 0.0;
else
{
f= u[i][i];
g= Math.sqrt(s);
if (f >= 0.0) g= -g;
h= f*g-s
u[i][i]=f-g;
for (j=l; j < n; j++)
{
s= 0.0
for (k=i; k < m; k++)
s += u[k][i]*u[k][j]
f= s/h
for (k=i; k < m; k++)
u[k][j]+=f*u[k][i]
}
}
q[i]= g
s= 0.0
for (j=l; j < n; j++)
s= s + u[i][j]*u[i][j]
if (s <= tolerance)
g= 0.0
else
{
f= u[i][i+1]
g= Math.sqrt(s)
if (f >= 0.0) g= -g
h= f*g - s
u[i][i+1] = f-g;
for (j=l; j < n; j++) e[j]= u[i][j]/h
for (j=l; j < m; j++)
{
s=0.0
for (k=l; k < n; k++)
s += (u[j][k]*u[i][k])
for (k=l; k < n; k++)
u[j][k]+=s*e[k]
}
}
y= Math.abs(q[i])+Math.abs(e[i])
if (y>x)
x=y
}
// accumulation of right hand gtransformations
for (i=n-1; i != -1; i+= -1)
{
if (g != 0.0)
{
h= g*u[i][i+1]
for (j=l; j < n; j++)
v[j][i]=u[i][j]/h
for (j=l; j < n; j++)
{
s=0.0
for (k=l; k < n; k++)
s += u[i][k]*v[k][j]
for (k=l; k < n; k++)
v[k][j]+=(s*v[k][i])
}
}
for (j=l; j < n; j++)
{
v[i][j] = 0;
v[j][i] = 0;
}
v[i][i] = 1;
g= e[i]
l= i
}
// accumulation of left hand transformations
for (i=n-1; i != -1; i+= -1)
{
l= i+1
g= q[i]
for (j=l; j < n; j++)
u[i][j] = 0;
if (g != 0.0)
{
h= u[i][i]*g
for (j=l; j < n; j++)
{
s=0.0
for (k=l; k < m; k++) s += u[k][i]*u[k][j];
f= s/h
for (k=i; k < m; k++) u[k][j]+=f*u[k][i];
}
for (j=i; j < m; j++) u[j][i] = u[j][i]/g;
}
else
for (j=i; j < m; j++) u[j][i] = 0;
u[i][i] += 1;
}
// diagonalization of the bidiagonal form
prec= prec*x
for (k=n-1; k != -1; k+= -1)
{
for (var iteration=0; iteration < itmax; iteration++)
{// test f splitting
var test_convergence = false
for (l=k; l != -1; l+= -1)
{
if (Math.abs(e[l]) <= prec){
test_convergence= true
break
}
if (Math.abs(q[l-1]) <= prec)
break
}
if (!test_convergence){
// cancellation of e[l] if l>0
c= 0.0
s= 1.0
var l1= l-1
for (i =l; i<k+1; i++)
{
f= s*e[i]
e[i]= c*e[i]
if (Math.abs(f) <= prec)
break
g= q[i]
h= pythag(f,g)
q[i]= h
c= g/h
s= -f/h
for (j=0; j < m; j++)
{
y= u[j][l1]
z= u[j][i]
u[j][l1] = y*c+(z*s)
u[j][i] = -y*s+(z*c)
}
}
}
// test f convergence
z= q[k]
if (l== k){
//convergence
if (z<0.0)
{ //q[k] is made non-negative
q[k]= -z
for (j=0; j < n; j++)
v[j][k] = -v[j][k]
}
break //break out of iteration loop and move on to next k value
}
console.assert(iteration < itmax-1, 'Error: no convergence.');
// shift from bottom 2x2 minor
x= q[l]
y= q[k-1]
g= e[k-1]
h= e[k]
f= ((y-z)*(y+z)+(g-h)*(g+h))/(2.0*h*y)
g= pythag(f,1.0)
if (f < 0.0)
f= ((x-z)*(x+z)+h*(y/(f-g)-h))/x
else
f= ((x-z)*(x+z)+h*(y/(f+g)-h))/x
// next QR transformation
c= 1.0
s= 1.0
for (i=l+1; i< k+1; i++)
{
g = e[i]
y = q[i]
h = s*g
g = c*g
z = pythag(f,h)
e[i-1] = z
c = f/z
s = h/z
f = x*c+g*s
g = -x*s+g*c
h = y*s
y = y*c
for (j =0; j < n; j++)
{
x = v[j][i-1]
z = v[j][i]
v[j][i-1] = x*c+z*s
v[j][i] = -x*s+z*c
}
z = pythag(f,h)
q[i-1] = z
c = f/z
s = h/z
f = c*g+s*y
x = -s*g+c*y
for (j =0; j < m; j++)
{
y = u[j][i-1]
z = u[j][i]
u[j][i-1] = y*c+z*s
u[j][i] = -y*s+z*c
}
}
e[l] = 0.0
e[k] = f
q[k] = x
}
}
// vt = transpose(v)
// return (u,q,vt)
for (i=0;i<q.length; i++)
if (q[i] < prec) q[i] = 0
// sort eigenvalues
for (i=0; i< n; i++){
// writeln(q)
for (j=i-1; j >= 0; j--){
if (q[j] < q[i]){
// writeln(i,'-',j)
c = q[j]
q[j] = q[i]
q[i] = c
for (k=0;k<u.length;k++) { temp = u[k][i]; u[k][i] = u[k][j]; u[k][j] = temp; }
for (k=0;k<v.length;k++) { temp = v[k][i]; v[k][i] = v[k][j]; v[k][j] = temp; }
i = j
}
}
}
return { U:u, S:q, V:v }
}
function pca(X,npc){
var USV = svd(X);
var U = USV.U;
var S = diag(USV.S);
var V = USV.V;
// T = X*V = U*S
var pcXV = dot(X,V)
var pcUdS = dot(U,S);
var prod = trunc(sub(pcXV,pcUdS), 1e-12);
var zero = zeros(prod.length, prod[0].length);
console.assert(same(prod,zero), 'svd and eig ways must be the same.');
return pcUdS;
}
};