Graph-tool 9

아래 예시 코드는 graph_tool 및 matplotlib 등을 사용하여 ers(블록 간 에지행렬), 노드 차수(node degrees), 그리고 fugacities 등을 시각화·분석하는 방법을 보여줍니다.
(사용 편의를 위해 단계별로 주석과 함께 제시하니, 필요한 부분만 골라서 쓰시면 됩니다.)

---

전체 코드 개요

1. 필요 라이브러리 임포트
2. 네트워크 불러오기 및 전처리: (가장 큰 연결성분, 자기루프·중복간선 제거)
3. 블록모델 추정 (minimize_blockmodel_dl)
4. ers(블록간 에지행렬) 확인 및 시각화
5. 노드 차수(node degrees) 분포·상세 정보 확인
6. fugacities(theta_out, mrs) 해석 및 시각화

---

import graph_tool.all as gt
import numpy as np
import matplotlib.pyplot as plt

###############################
# 1) 네트워크 로드 및 전처리
###############################

# (예시) 정치 블로그 데이터 로드
g_orig = gt.collection.data["polblogs"]

# (a) 가장 큰 연결성분으로 필터링
g_largest = gt.GraphView(g_orig, vfilt=gt.label_largest_component(g_orig))

# (b) 복제 & 단순그래프 만들기 (자기루프/중복간선 제거)
g = gt.Graph(g_largest, prune=True)
gt.remove_self_loops(g)
gt.remove_parallel_edges(g)

print("== After filtering & removing loops ==")
print("Num of vertices:", g.num_vertices())
print("Num of edges   :", g.num_edges())

###############################
# 2) 블록 모델 추정
###############################

# 블록모델 추정 (MDL 기반)
state = gt.minimize_blockmodel_dl(g)

# 블록(그룹) 식별자 (각 노드가 어느 블록에 속하는지)
b_array = state.b.a

# 블록 개수 확인
unique_blocks = np.unique(b_array)
print("== Block model info ==")
print("Number of blocks:", len(unique_blocks))

###############################
# 3) ers(블록 간 에지행렬) 확인
###############################

# state.get_bg() : 블록 단위 메타그래프
# state.get_ers(): 블록 간 엣지 수(혹은 가중치)
ers = gt.adjacency(state.get_bg(), state.get_ers()).T
# ers.shape -> (B, B) 형태 (B = 블록 수)

print("== ers (block-block edge matrix) shape:", ers.shape)

# (a) 대각합, 전체합
diag_sum = np.trace(ers)        # 대각(블록 내부) 합
off_diag_sum = np.sum(ers) - diag_sum
print(f"Sum of diagonal (intra-block): {diag_sum}")
print(f"Sum of off-diagonal (inter-block): {off_diag_sum}")

# (b) 만약 시각화를 하고 싶다면:
plt.figure(figsize=(6,5))
plt.imshow(ers, cmap='Blues', origin='lower')
plt.colorbar(label='Number of edges')
plt.title("Block-Block Edge Matrix (ers)")
plt.xlabel("Block ID")
plt.ylabel("Block ID")
plt.tight_layout()
plt.show()

# --> 이 행렬을 보고, 대각선 성분이 큰지, 어떤 특정 블록쌍이 유난히 큰지 확인 가능
# --> ex) 대각선이 압도적으로 크다면 "assortative-like" / 비대각이 두드러지면 "bipartite/core-periphery" 가능성.

###############################
# 4) 노드 차수(node degrees) 분석
###############################

# g가 무향 그래프라 가정 시
degrees = g.degree_property_map("out").a  # "in"도 같음

print("== Node degree info ==")
print("Min degree:", np.min(degrees))
print("Max degree:", np.max(degrees))
print("Average degree:", np.mean(degrees))

# 차수 분포 히스토그램
plt.figure(figsize=(6,4))
plt.hist(degrees, bins=30, edgecolor='black', alpha=0.7)
plt.xlabel("Degree")
plt.ylabel("Count of nodes")
plt.title("Degree distribution")
plt.tight_layout()
plt.show()

# (a) 만약 누가 허브인지 궁금하다면:
top_k = 5
sorted_idx = np.argsort(-degrees)  # 내림차순 정렬 (큰 차수부터)
print(f"Top {top_k} highest-degree nodes:")
for rank in range(top_k):
    node_id = sorted_idx[rank]
    print(f"  Rank {rank+1} node={node_id}, degree={degrees[node_id]}")

###############################
# 5) Fugacities (theta_out, mrs 등) 계산
###############################

# (a) 블록간 실제 에지수(ers), 노드별 out/in degree 가져오기
ers_matrix = ers  # 이미 구함
out_degs = g.degree_property_map("out").a
in_degs  = g.degree_property_map("in").a  # 무향인 경우 out_degs와 동일

# (b) solve_sbm_fugacities 호출
mrs, theta_out, theta_in = gt.solve_sbm_fugacities(b_array, ers_matrix, out_degs, in_degs)

print("== Fugacities (theta) & mrs matrix ==")
print("  mrs shape:", mrs.shape)  # (B, B)
print("  theta_out shape:", theta_out.shape)  # (N,)
# 만약 무향이면 theta_in == theta_out 비슷할 것

# (c) mrs: 블록 r-s 간 fugacity
# 대략적으로 "r, s 블록 사이 연결 발생 세기"라 보면 됨
# (참고: ers와 mrs는 비례관계지만, 정규화/기대값 설정이 조금 다름)

# 만약 시각화 해보고 싶다면:
plt.figure(figsize=(6,5))
plt.imshow(mrs, cmap='Reds', origin='lower')
plt.colorbar(label='Fugacity (mrs[r,s])')
plt.title("Block-Block fugacity matrix (mrs)")
plt.xlabel("Block ID")
plt.ylabel("Block ID")
plt.tight_layout()
plt.show()

# (d) theta_out: 각 노드별 fugacity
# 이 값이 클수록, 모델에서 "연결을 많이 유발"하는 노드
top_t = 5
sorted_theta_idx = np.argsort(-theta_out)
print(f"Top {top_t} largest theta_out nodes:")
for rank in range(top_t):
    node_id = sorted_theta_idx[rank]
    print(f"  Rank {rank+1} node={node_id}, theta_out={theta_out[node_id]:.4f}, degree={degrees[node_id]}")

# --> 실제 차수와 theta_out 간 비교해보면
#     "theta_out 큰 노드 = 차수 큰 노드" 패턴이 보일 가능성 높음.

###############################
# (Optional) 추가 해석
###############################
# - ers와 mrs 모두 (B,B) 매트릭스이므로, 이를 종합해:
#    1) ers 대각선이 큰지/off-diagonal이 큰지
#    2) mrs 대각선이 큰지/off-diagonal이 큰지
#      => 어떤 블록들이 자기 내부 연결이 세거나, 어느 블록 쌍이 강력히 연결됐는지
# - theta_out 상위 노드들이, 실제로 네트워크 전체 허브인지 체크.
# - 차수 분포(히스토그램)에서, 파워법 모양인지(직선로그-로그?), 정규형인지, etc. 
# - 이런 분석을 통해 네트워크 구조 해석 가능!

---

코드 해석 및 활용 팁

1. ers(블록 간 에지행렬)의 시각화

   • plt.imshow(ers, ...) 부분에서 파란(또는 다른 컬러맵) 색깔이 진한 블록-블록 쌍은 그 둘 사이에 에지가 많다는 뜻.
   • 대각선이 유독 진하면 블록 내부 연결이 많으니, 전형적인 “assortative(내부 결속)” 구조를 의심.
   • 특정 비대각 구간만 진하면 “bipartite / core-periphery” 등 다른 구조 가능성.

2. 노드 차수(degrees) 히스토그램

   • 고차수 노드가 몇 개만 튀어나온다면, 허브가 존재한다고 해석.
   • 차수분포가 치우쳐 있으면(예: 파워법) “hub-and-spoke” 형태일 수도 있고, 균일하면 “random-like”일 수도.

3. fugacities

   • theta_out이 큰 노드는 모델에서 “많은 연결을 일으키는(=고차수)” 노드로 본다는 의미.
   • mrs[r, s]가 큰 블록 쌍은 “연결 강도가 크다” → 해당 블록 쌍에서 에지가 자주 생김.

4. 추가로

   • 커뮤니티 라벨 b_array를 활용해, “같은 블록 노드들만 추출”하여 서브그래프 시각화 등도 가능.
   • 더 높은 수준 분석(ex. nested SBM, modular hierarchy)으로 발전할 수도 있음.

---

위 코드를 실행하면,

• 블록 간 에지행렬(ers) 시각화,
• 노드 차수 분포,
• fugacity(mrs, theta_out) 값
  등을 구체적으로 확인하여, 네트워크 구조(커뮤니티·허브·assortative 등)를 좀 더 깊이 해석할 수 있게 됩니다.

아래 튜토리얼은 (\texttt{graph_tool.generation.solve_sbm_fugacities}) 함수가 무엇을 하는지, 그리고 이를 어떻게 활용하여 SBM(확률적 블록 모델)의 'fugacity'(퓨가시티, 일종의 라그랑주 승수와 유사한 개념)를 효율적으로 추정할 수 있는지 알려주는 “초보자를 위한 가이드”입니다.

• 핵심 목표: 블록(그룹) 간 예상 엣지 수(ers)와 각 노드의 기대 차수(out_degs, in_degs)가 주어졌을 때, ( \text{fugacity} )들을 계산하는 방법을 이해하고, 실제 코드를 통해 다뤄보기.

---

1. 이 함수가 해결해주는 문제는?

(1) SBM에서 “fugacity”란?

• SBM(확률적 블록 모델)은 네트워크를 여러 그룹(블록)으로 나누고, 그룹-그룹 간 연결 확률(또는 개수)을 통해 네트워크를 생성하는 모델입니다.
• “fugacity”(퓨가시티)는 수학적으로 라그랑주 승수(혹은 지수화한 형태)로 해석할 수 있으며,
  • 노드 차수(degree) 쪽 fugacity: 각 노드가 몇 개의 엣지를 기대하는지 조절한다.
  • 그룹 간 edge count 쪽 fugacity: 어떤 두 블록 사이에 몇 개의 엣지가 형성되는지를 조절한다(단순그래프 / 멀티그래프 구분에 따라 세부 식이 달라짐).

(2) solve_sbm_fugacities(...) 함수의 의도

• 입력:
  1. b: 각 노드가 어느 블록에 속하는지 (길이 (N)의 배열)
  2. ers: 블록-블록 간 “예상(평균) 엣지 수” 행렬 ((B\times B))
  3. out_degs / in_degs: 각 노드의 “예상 out-degree” / “예상 in-degree”
     • 무향 그래프일 경우 in_degs 없이도 됨(동일하다고 간주).
  4. multigraph, self_loops: 멀티그래프인지, 자기 루프 허용인지 등.
• 출력:
  1. mrs: 블록-블록 간 fugacity 행렬 ((B\times B)),
  2. out_theta, (옵션) in_theta: 각 노드별 fugacity(길이 (N)).
• 이들은 SBM에서 “최대 엔트로피 조건 + 주어진 기대값”이 성립하도록 자기일관성 방정식을 해석학적으로(또는 수치 해법으로) 풀어서 얻어지는 해.

(3) 어떤 경우에 사용하는가?

1. 이미 블록분할(e.g. b_array)이 정해져 있고, 블록-블록 간 엣지 개수 기대치((ers))와 노드별 차수(기대치)가 주어져 있을 때, 최대 엔트로피 SBM(또는 Poisson SBM)에 필요한 파라미터((\lambda_{rs}) 등)를 추정하고 싶다.
2. 이후 이 fugacity들을 기반으로 (\texttt{generate_maxent_sbm}) 함수 등을 써서 랜덤 그래프 샘플링을 할 수 있음.
   • 즉, “실제 네트워크와 동일한 블록 구조와 노드 차수를 가진 (혹은 비슷한) 최대 엔트로피 그래프를 생성”할 수 있다.

---

2. 이론적으로는?

문서에 따르면, 단순그래프(undirected)일 때나 유향/멀티그래프일 때, 수학적 자기일관성 식은 각각 다릅니다.

• 예를 들어, 무향 단순그래프에서 각 블록 쌍 ((r,s))에 대해
  [
  \mathbb{E}[A_{ij}] = \frac{ \thetai \theta_j m{rs} }{ 1 + \thetai \theta_j m{rs} },
  ]
  등등이 성립해야 하고, 여기서
  • (m_{rs})가 “그룹 (r)와 (s) 사이의 (\text{fugacity})”
  • (\theta_i)가 “노드 (i)의 (\text{fugacity})”.

이 함수는 내부적으로 “이 방정식을 만족”하도록 (\thetai)와 (m{rs})를 찾아주는 과정을 진행합니다.

---

3. 주요 파라미터 설명

문서에서 발췌 & 정리:

1. b:
   • 길이 (N)의 정수 배열. b[i]는 노드 (i)가 어느 블록(r)에 속하는지.
2. ers:
   • 블록 간 expected edge count 행렬, 크기 ((B,B)).
   • ers[r,s] = 그룹 (r)-(s) 사이에 평균 몇 개의 엣지(또는 (r=s)면 내부엣지 2배 계산)
3. out_degs, in_degs:
   • 노드별 기대 차수(유향이면 in/out 다를 수 있음).
   • 무향이면 in_degs=None으로 둬도 됨.
4. multigraph, self_loops:
   • 병렬엣지(멀티엣지) 허용 여부, 자기 루프 허용 여부.
   • 수학적 방정식에 영향을 줌.
5. epsilon: 수치해석 상 “정밀도” 임계값.
6. iter_solve: True면 단순 반복법, False면 scipy의 minimize/root solver 쓸 수도 있음(하지만 보통 True).
7. max_iter: 반복 횟수 제한(0이면 무제한).
8. verbose: 진행 상황 로깅.

---

4. 간단 예제 코드 (유사 스니펫)

아래는 “가정된 예시 데이터”를 이용해 튜토리얼 형태로 작성한 코드 예시입니다.
주요 흐름:
1. 임의의 블록 분할 (b), 블록 엣지행렬 ers, 노드 차수 out_degs 등을 준비
2. solve_sbm_fugacities(...)로 (mrs, theta_out, theta_in) 계산
3. 결과를 해석

import numpy as np
import graph_tool.all as gt

def demo_solve_sbm_fugacities():
    """
    demo_solve_sbm_fugacities:
      - 임의의 블록 구조, ers, 노드 차수 등을 준비해, 
        solve_sbm_fugacities 함수를 호출해본다.
    """
    
    # 1) 블록 분할 (예: 노드가 6개 있다고 가정)
    b = np.array([0, 0, 1, 1, 2, 2], dtype=int)  
    # => 6노드를 3개 블록(0,1,2)로 나눈 예시
    
    # 2) ers(블록 간 에지 기대치) 준비
    #    B=3 블록, ers[r,s] = 그룹r-s의 예상 엣지 수
    #    예: 대각성분은 블록 내부엣지
    ers_data = np.array([
        [10.0,  3.0,  1.0],
        [ 3.0, 12.0,  2.0],
        [ 1.0,  2.0,  5.0]
    ])
    # => r=0 블록 내부 10개, r=1 블록 내부 12개 등
    
    # 3) out_degs (노드 차수) 준비
    #    여기서는 6노드 각각 임의로 설정
    out_degs = np.array([5, 4, 6, 6, 2, 3], dtype=float)
    # in_degs: 무향 그래프라면 None으로도 가능
    in_degs = None  
    
    # 4) solve_sbm_fugacities 호출
    mrs, theta_out = gt.solve_sbm_fugacities(
        b, ers_data, out_degs, in_degs, 
        multigraph=False,   # 단순그래프라고 가정
        self_loops=False,   # 자기루프X
        epsilon=1e-8,
        iter_solve=True,    # 단순 반복법
        verbose=True
    )
    # => in_degs=None 이므로 반환값은 (mrs, theta_out)
    #    in_degs != None 이면 (mrs, theta_out, theta_in) 세 개를 반환
    
    print("[결과] mrs shape:", mrs.shape)   # (3,3)
    print("[결과] theta_out shape:", theta_out.shape) # (6,)
    
    # 5) 결과 해석
    # mrs: (B,B) sparse matrix. 필요하다면 .toarray()로
    mrs_dense = mrs.toarray()
    print("[결과] mrs (dense):\n", mrs_dense)
    
    # theta_out: 노드별 퓨가시티 => 높을수록 엣지를 많이 생성하려는 경향
    for i, val in enumerate(theta_out):
        print(f"  Node {i}, theta_out={val:.4f}, expectedDeg={out_degs[i]}")

if __name__ == "__main__":
    demo_solve_sbm_fugacities()

(실행 결과 예시)

[결과] mrs shape: (3, 3)
[결과] theta_out shape: (6,)
[결과] mrs (dense):
 [[ 1.2839  0.7542  0.1031]
  [ 0.7542  1.4943  0.2205]
  [ 0.1031  0.2205  0.7127]]
  Node 0, theta_out=0.5332, expectedDeg=5.0
  Node 1, theta_out=0.4928, expectedDeg=4.0
  Node 2, theta_out=0.6211, expectedDeg=6.0
  Node 3, theta_out=0.6026, expectedDeg=6.0
  Node 4, theta_out=0.2113, expectedDeg=2.0
  Node 5, theta_out=0.3201, expectedDeg=3.0

(단, 위 수치는 가상의 예시로서, 실제와 다름.)

---

5. 결과물 해석 포인트

1. mrs[r, s](블록 간 fugacity):
   • 값이 클수록 블록 (r)-(s) 간 연결이 더 활발히 일어날 가능성 크다.
   • (B,B) 대칭행렬(무향) 혹은 비대칭(유향)일 수 있음.
2. theta_out[i](노드 i의 fugacity):
   • 노드 i가 “얼마나 많은 엣지를 만들어내는가”에 대한 기여도.
   • 차수가 큰 노드는 보통 theta_out[i]도 큼.
3. 무향 + 자기루프 X 인 경우, 식이 “(\frac{\thetai\theta_j \cdot m{rs}}{1 + \thetai\theta_j m{rs}})” 형태” 등으로 포화가 발생.
   • 실제로 solve_sbm_fugacities가 이 식을 만족하도록 (\theta, m)들을 찾는다.

---

6. 주의 사항 / 자주 겪는 문제

1. 블록 분할(b)가 잘못됨:
   • b는 “노드 i가 속한 블록”을 나타내는 배열. 즉, 0부터 (B-1)까지 정수 값을 가져야 하며, 블록 수가 몇 개인지와 일치해야 한다.
2. ers의 shape:
   • (B,B)여야 함. 가끔 (N,N)짜리 adjacency와 혼동해 잘못 넘기는 경우가 있으니 주의.
3. in_degs가 None인지 아닌지에 따라 반환값이 달라짐. (유향 vs 무향)
4. 수렴 문제:
   • 매우 큰 네트워크이거나, 특정 노드 차수가 엄청 크면, 방정식이 잘 안맞아 수렴 어려울 수 있음 => epsilon, max_iter 조절.
5. multigraph, self_loops 옵션:
   • True로 하면 “Poisson-like” 형태. False면 simple graph 형태로 인한 “포화(saturation)”가 발생.

---

7. 결론

• (\texttt{solve_sbm_fugacities})는 SBM(최대 엔트로피 or Poisson) 그래프를 생성할 때 필요한 라그랑주 승수(=fugacity)를 자동으로 계산해주는 강력한 도구입니다.
• 사전에 블록 분할과 원하는 “평균 블록 간 엣지 수(ers)” 그리고 “노드별 차수(예상)”를 지정해두면, 이 함수가 (\theta)와 (m_{rs})를 맞춰줍니다.
• 이후 generate_maxent_sbm 함수를 통해 샘플 그래프를 만들거나, 혹은 모델 해석(고차수 노드, 특정 블록들 연결성) 등 다양한 분석을 할 수 있습니다.

추천 후속 학습

• (\texttt{generate_maxent_sbm})로 실제 네트워크 표본을 생성하는 코드 실습
• 다양한 파라미터(멀티그래프 허용 여부, 자기루프 등)에 따른 fugacity 값 변화를 관찰
• 실제 네트워크 데이터를 가지고, 블록모델 추정 → solve_sbm_fugacities 수행 → “원본과 비슷한 확률모형” 만들어 보기.