import numpy as np
import scipy
from scipy import sparse
import array
#Brute force function to check for all 3 cliques by checking
#mutual neighbors between 3 vertices
def get3CliquesBrute(Edges):
[I, J, V] = [[], [], []]
NVertices = len(Edges) #返回对象（字符、列表、元组等）长度或项目个数
MaxNum = int(NVertices*(NVertices-1)*(NVertices-2)/6)
J = np.zeros((MaxNum, 3))
[i, j, k] = [0, 0, 0]
edgeNum = 0
for i in range(NVertices):#可创建一个整数列表
for j in Edges[i]:
if j < i:
continue #跳出本次循环
for k in Edges[j]:
if k < j or k < i:
continue
if k in Edges[i]:
[a, b, c] = sorted([i, j, k])
J[edgeNum, :] = [Edges[a][b], Edges[a][c], Edges[b][c]]
edgeNum += 1
J = J[0:edgeNum, :]
V = np.zeros(J.shape)
V[:, 0] = 1
V[:, 1] = -1
V[:, 2] = 1
I = np.zeros(J.shape)
for k in range(3):
I[:, k] = np.arange(I.shape[0])
return (I, J, V)
#Extract 3 cliques from maximal cliques, given an array
#with edge indices
def get3CliquesFromMaxCliques(Cliques, E):
cliques = set()#建一个无序不重复元素集，可进行关系测试，删除重复数据，还可以计算交集、差集、并集等
[I, J, V] = [[], [], []]
for c in Cliques:
for i in range(len(c)):
for j in range(i+1, len(c)):
for k in range(j+1, len(c)):
cliques.add((c[i], c[j], c[k]))
for c in cliques:
I.append(3*[len(I)])
J.append([E[c[0], c[1]]-1, E[c[0], c[2]]-1, E[c[1], c[2]]-1])
V.append([1.0, -1.0, 1.0])
return (I, J, V)
#Compare boundary matrices up to a permutation of the rows
def compareBoundaryMatricesModPerm(A, B):
setA = set()
for i in range(A.shape[0]): #Assuming csr
a = A[i, :]
a = a.tocoo()
arr = tuple(a.col.tolist() + a.data.tolist())
setA.add(arr)
setB = set()
for i in range(B.shape[0]):
b = B[i, :]
b = b.tocoo()
arr = tuple(b.col.tolist() + b.data.tolist())
setB.add(arr)
if len(setA.difference(setB)) > 0 or len(setB.difference(setA)) > 0:
return False
return True
```python
import numpy as np
import matplotlib.pyplot as plt
import scipy.io as sio
import scipy
from scipy import sparse
from scipy.sparse.linalg import lsqr, cg, eigsh
import time
import jellyfish
from CliqueAlgorithms import *
def getKendallTau(order1, order2):
"""
Given two global rankings, return the Kendall Tau Score
"""
N = len(order1)
rank1 = np.zeros(N)
rank1[order1] = np.arange(N)
rank2 = np.zeros(N)
rank2[order2] = np.arange(N)
A = np.sign(rank1[None, :] - rank1[:, None])
B = np.sign(rank2[None, :] - rank2[:, None])
return np.sum(A*B)/float(N*(N-1))
#tau, p_value = scipy.stats.kendalltau(rank1, rank2)
def getJWDistance(order1, order2):
"""
Given two global rankings, return the Jaro Winkler Distance
"""
s = u"abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ"
s1 = u""
s2 = u""
for i in range(len(order1)):
s1 += s[order1[i]]
for i in range(len(order2)):
s2 += s[order2[i]]
return jellyfish.jaro_winkler(s1, s2)
def makeDelta0(R):
"""
Return the delta0 coboundary matrix
:param R: NEdges x 2 matrix specifying edges, where orientation
is taken from the first column to the second column
R specifies the "natural orientation" of the edges, with the
understanding that the ranking will be specified later
It is assumed that there is at least one edge incident
on every vertex
"""
NVertices = int(np.max(R)) + 1
NEdges = R.shape[0]
#Two entries per edge
I = np.zeros((NEdges, 2))
I[:, 0] = np.arange(NEdges)
I[:, 1] = np.arange(NEdges)
I = I.flatten()
J = R[:, 0:2].flatten()
V = np.zeros((NEdges, 2))
V[:, 0] = -1
V[:, 1] = 1
V = V.flatten()
Delta = sparse.coo_matrix((V, (I, J)), shape=(NEdges, NVertices)).tocsr()
return Delta
def makeDelta1(R, verbose):
"""Make the delta1 coboundary matrix
:param R: Edge list NEdges x 2. It is assumed that
there is at least one edge incident on every vertex
"""
NEdges = R.shape[0]
NVertices = int(np.max(R))+1
#Make a list of edges for fast lookup
Edges = []
for i in range(NVertices):
Edges.append({})
for i in range(R.shape[0]):
[a, b] = [int(R[i, 0]), int(R[i, 1])]
Edges[a][b] = i
Edges[b][a] = i
tic = time.time()
(I, J, V) = get3CliquesBrute(Edges)
toc = time.time()
if verbose:
print("Elapsed time 3 cliques brute: %g"%(toc - tic))
[I, J, V] = [a.flatten() for a in [I, J, V]]
TriNum = int(len(I)/3)
Delta1 = sparse.coo_matrix((V, (I, J)), shape = (TriNum, NEdges)).tocsr()
return Delta1
def doHodge(R, W, Y, verbose = False):
"""
Given
:param R: NEdges x 2 matrix specfiying comparisons that have been made
:param W: A flat array of NEdges weights parallel to the rows of R
:param Y: A flat array of NEdges specifying preferences
:returns: (s, I, H): s is scalar function, I is local inconsistency vector,
H is global inconsistency vector
"""
#Step 1: Get s
if verbose:
print("Making Delta0...")
tic = time.time()
D0 = makeDelta0(R)
toc = time.time()
if verbose:
print("Elapsed Time: %g seconds"%(toc-tic))
wSqrt = np.sqrt(W).flatten()
WSqrt = scipy.sparse.spdiags(wSqrt, 0, len(W), len(W))
WSqrtRecip = scipy.sparse.spdiags(1/wSqrt, 0, len(W), len(W))
A = WSqrt*D0
b = WSqrt.dot(Y)
s = lsqr(A, b)[0]
#Step 2: Get local inconsistencies
if verbose:
print("Making Delta1...")
tic = time.time()
D1 = makeDelta1(R, verbose)
toc = time.time()
if verbose:
print("Elapsed Time: %g seconds"%(toc-tic))
B = WSqrtRecip*D1.T
resid = Y - D0.dot(s) #This has been verified to be orthogonal under <resid, D0*s>_W
u = wSqrt*resid
if verbose:
print("Solving for Phi...")
tic = time.time()
Phi = lsqr(B, u)[0]
toc = time.time()
if verbose:
print("Elapsed Time: %g seconds"%(toc-tic))
I = WSqrtRecip.dot(B.dot(Phi)) #Delta1* dot Phi, since Delta1* = (1/W) Delta1^T
#Step 3: Get harmonic cocycle
H = resid - I
return (s, I, H)
def getWNorm(X, W):
return np.sqrt(np.sum(W*X*X))
```python
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from Hodge import *
def getConsistencyRatios(Y, I, H, W, verbose=False):
normD0s = getWNorm(Y-H-I, W)
[normY, normI, normH] = [getWNorm(Y, W), getWNorm(I, W), getWNorm(H, W)]
a = (normD0s/normY)**2
b = (normI/normY)**2
c = (normH/normY)**2
if verbose:
print("|D0s/Y| = %g"%a)
print("Local Inconsistency = %g"%b)
print("Global Inconsistency = %g"%c)
print("a + b + c = %g"%(a + b + c))
return (a, b, c)
def doTotalOrderExperiment(N, NComparisons, NQuant = -1):
"""
Do an experiment comparing binary weights on a total order
to real weights
Parameters
----------
N: int
Number of items to consider
NComparisons: int
Number of comparisons performed
NQuant: int
Number of quantization levels
"""
I, J = np.meshgrid(np.arange(N), np.arange(N))
I = I[np.triu_indices(N, 1)]
J = J[np.triu_indices(N, 1)]
#[I, J] = [I[0:N-1], J[0:N-1]]
NEdges = len(I)
R = np.zeros((NEdges, 2))
R[:, 0] = J
R[:, 1] = I
#W = np.random.rand(NEdges)
W = np.ones(NEdges)
#Note: When using binary weights, Y is not necessarily a cocycle
Y = I - J
if NQuant > -1:
Y = NQuant*np.round(NQuant*Y/N)
idx = np.random.permutation(Y.size)[0:NComparisons]
R = R[idx, :]
W = W[idx]
Y = Y[idx]
(s, I, H) = doHodge(R, W, Y, verbose=False)
return getConsistencyRatios(Y, I, H, W)
if __name__ == '__main__2':
np.random.seed(0)
N = 50
NTrials = 10
AllNComparisons = np.arange(100, 1000, 10)
Is = np.zeros((AllNComparisons.size, NTrials))
Hs = np.zeros_like(Is)
NLevels = 3
for i, NComparisons in enumerate(AllNComparisons):
for j in range(NTrials):
d0, IMag, HMag = doTotalOrderExperiment(N, NComparisons, NLevels)
Is[i, j] = IMag
Hs[i, j] = HMag
print("%i, local:%.3g, global:%.3g"%(NComparisons, np.mean(Is[i, :]), np.mean(Hs[i, :])))
AllNComparisons = AllNComparisons[:, None]*np.ones((1, NTrials))
sns.lineplot(AllNComparisons.flatten(), Is.flatten())
sns.lineplot(AllNComparisons.flatten(), Hs.flatten())
plt.xlabel("Number of Comparisons")
plt.ylabel("Inconsistency Ratio")
plt.legend(["Local Inconsistencies", "Global Inconsistencies"])
plt.title("%i Level Inconsistencies for %i Objects"%(NLevels, N))
plt.show()
if __name__ == '__main__':
import glob#glob模块提供了一个函数用于从目录通配符搜索中生成文件列表
files = glob.glob("Results/*")#匹配所有的符合条件的文件，并将其以list的形式返回。
for f in files:
#R = sio.loadmat('R.mat')['R']
R = np.loadtxt(f)
print(f, R.shape[0])
[R, Y] = [R[:, 0:2], R[:, 2]]
W = np.random.rand(len(Y))
#W = np.ones(len(Y))
(s, I, H) = doHodge(R, W, Y)
print(np.argsort(s))
getConsistencyRatios(Y, I, H, W, verbose=True)
