요통, 엉덩이 통증, 무릎 통증, 보행 이상 간의 양방향 관계 조사: 멘델 무작위화 연구

요약

IVW 분석 결과, 허리 통증과 엉덩이 통증이 보행 이상과 연관될 수 있으나, 무릎 통증과 역연관성은 확실하지 않습니다.

초록

보행 이상은 하체 생체역학 변화와 근골격계 질환에서 기능적 장애와 관련이 있습니다. 관찰 연구에서 요통, 엉덩이 통증, 무릎 통증이 보행 장애와 연관이 있음에도 불구하고, 이들의 인과관계는 명확하지 않습니다. 우리는 부위 특이적 통증과 보행 이상 간의 잠재적 인과관계를 평가하기 위해 양방향 멘델 무작위화(MR) 연구를 수행하였습니다. 대규모 유전체 연관 연구에서 독립적인 단일 염기체 다형성을 요통, 엉덩이 통증, 무릎 통증의 유전적 도구로 사용했습니다. 보행 이상은 출처 GWAS에서 자기보고된 걷기 어려움으로 정의되었으며, 전방 MR 분석과 역방향 MR 분석에서 노출 결과로 사용되었습니다. 역분산 가중(IVW) 방법이 주요 분석이었으며, 이질성, 수평 다면성, 강인성에 대한 민감도 분석이 이를 뒷받침하였다. 전방 분석에서 유전적으로 예측된 요통(OR = 1.53, 95% CI = 1.058–2.219, P = 0.024)과 고관절 통증(OR = 1.55, 95% CI = 1.067–2.252, P = 0.021)은 IVW 분석에서 보행 이상 위험 증가와 연관되어 있었습니다. 그러나 네 가지 추가 MR 방법은 통계적으로 유의미하지 않았으나, 추정치는 대체로 방향성 면에서 일관되었습니다. 유전적으로 예측된 무릎 통증은 유의미한 연관성이 없었으며(OR = 2.15, 95% 신뢰구간 = 0.234–19.789, P = 0.498), 넓은 신뢰 구간에서 상당한 불확실성이 반영되었습니다. 역방향 MR 분석에서는 보행 이상에 대한 유전적 위험성이 허리, 엉덩이, 무릎 통증 위험을 증가시킨다는 증거를 찾지 못했으나, 보행 이상 SNP는 4개만 발견되었습니다. 유의미한 이질성이나 수평 다면형성은 검출되지 않았습니다. 전반적으로 본 연구는 요통과 고관절통이 보행 이상과 연관될 수 있다는 제한된 IVW 기반 유전적 증거를 제공합니다. 발견은 신중하게 해석되어야 하며, 더 정제된 표현형을 가진 더 큰 GWAS 데이터셋을 사용해 검증해야 합니다.

서론

보행 이상은 운동 또는 감각 장애에서 비롯될 수 있으며, 임상 특징은 기저 병리의 위치와 성격에 따라 달라집니다. 특정한 비정상적인 보행 패턴은 특정 질병에 대한 중요한 진단 단서를 제공할 수 있습니다. 특히 근골격계 질환에서 보행 이상은 독특한 생체역학적 특성을 동반하며 통증, 관절 기능 장애, 운동성 저하와 밀접한 관련이 있습니다 ^1,2,3. 연구에 따르면 장기적인 보행 이상이 골관절염 환자의 상태를 악화시킬 수^{있습니다.} 더불어, 지속적인 비정상적인 보행은 운동 관절의 불안정성을 증가시키고 관절 부하를 변화시킬 수 있습니다; 비정상적인 역학 하에서 반복적인 관절 활동은 관절 손상을 더욱 악화시켜 통증을 유발하고 점진적인 근육 위축을 촉진하며 근력 저하를 더욱 심화시킬 수 있습니다. 이 요인들이 상호작용하여 악순환을 형성하여 근골격계 장애 환자에서 통증과 기능 장애가 점진적으로 악화되는 데 기여할 수 있습니다⁵. 따라서 부위 특이적 근골격계 통증과 보행 이상 간의 연관성 방향과 잠재적 인과관계를 명확히 하는 것은 관찰 통증-보행 관계를 해석하고 향후 보행 기능 장애 및 재활 연구를 이해하는 데 도움이 될 수 있습니다.

그러나 요통, 엉덩이 통증, 무릎 통증과 같은 부위 특이적 통증이 보행 수행 장애와 연관된 관찰 증거가 축적되었음에도 불구하고, 기존 관찰 연구는 교란, 측정 오차, 역인과관계에 취약하여 이러한 연관성의 방향과 인과관계를 규명하기 어렵습니다. 또한, 통증과 보행의 관계는 말초에서 수신되는 통각수용 정보의 중추신경계 적응과 조절을 포함할 수 있어, 관찰된 통증-보행 연관성 해석을 더욱 복잡하게 만듭니다. 멘델 무작위화(MR)는 관찰 연구의 주요 한계를 극복하기 위해 고안된 역학적 접근법으로, 인과 추론 연구에 널리 적용되어 왔습니다. MR은 독립적인 단일 염기화 다형성(SNP)을 도구 변수로 사용하여 노출과 결과 간의 잠재적 인과관계를 추론합니다^. 예를 들어, 개인에게 요통을 유발하도록 지정하는 대신, MR은 대규모 유전체 연관 연구(GWAS) 데이터셋에서 요통에 대한 유전적 취약성이 보행 관련 결과와 연관되어 있는지 조사할 수 있습니다. 유전적 변이가 수정 시점, 즉 배아 형성 시점에 무작위로 할당되기 때문에, MR은 나이, 체질량지수, 신체 활동, 그리고 동반되는 근골격계 질환과 같은 교란 요인에 의한 편향을 줄일 수 있습니다. 이 접근법은 강력한 GWAS 데이터와 유효한 유전 도구가 있을 때 집단 수준의 인과관계를 탐구하는 데 적합합니다. 그러나 본 연구에서 사용된 표현형은 광범위한 자기보고된 GWAS 특성으로, 통증의 심각도, 증상 지속 기간, 임상의가 확정한 진단, 객관적인 보행 매개변수를 완전히 포착하지 못할 수 있습니다. 수정 시 유전적 변이의 무작위 배분을 활용함으로써, MR은 교란 요인에 의한 편향을 효과적으로 완화하고 인과관계를 역전시켜 더 그럴듯한 인과관계를 도출할 수 있습니다.

본 연구의 목적은 GIWAS의 자가 보고된 요통, 엉덩이 통증, 무릎 통증에 대한 요약 데이터와 IEU 오픈 GWAS 데이터베이스의 보행 이상에 관한 GWAS 데이터를 사용하여 양방향 MR 분석을 수행하는 것이었습니다. 본 연구에서 통증 표현형은 광범위한 자기보고된 부위 특이적 통증 특성을 나타내며, 보행 이상은 임상의가 확인한 진단, 실험실에서 측정한 보행 지표, 표준화된 기능 점수보다는 자기보고 걷기 어려움을 나타냈습니다. 이 양방향 설계는 GWAS 요약 통계를 사용하여 부위 특이적 근골격계 통증과 보행 이상 간의 잠재적 인과관계 방향을 평가하는 것을 목표로^{했습니다.}

프로토콜

이 연구는 IEU OpenGWAS 데이터베이스와 GWAS 카탈로그에서 공개된 전유전체 연관성 연구(GWAS) 요약 수준 데이터를 사용했습니다. 개별 참가자 데이터는 접근하지 않았으며, 새로운 인간 참가자도 모집되지 않았습니다. 원래의 기여 GWAS 연구에서는 윤리적 승인과 서면 동의가 이루어졌습니다. 따라서 공개된 요약 통계의 2차 분석에 대해 추가적인 기관 검토 위원회의 승인은 필요하지 않았습니다.

연구 원칙과 설계

MR 연구의 기본 단계에는 GWAS 요약 데이터 획득, SNP 선별 및 평가, 통계 분석 수행, 품질 관리 조치 시행이 포함됩니다. 그림 1에서 보듯, MR 분석의 정확도는 세 가지 중요한 가정의 충족에 기반합니다: (1) 관련성 가정: IV는 노출 표현형, 즉 요통, 엉덩이 통증, 무릎 통증과 관련되어야 합니다⁸; (2) 독립성 가정: IV는 "노출-결과"에 영향을 미치는 교란 요인과 무관하며; (3) 배타성 가정: IV는 노출을 통해서만 결과에 영향을 미치고 다른 경로는 영향을 받지 않는다⁹. 노출 인자와 결과 변수는 각 분석에서 서로 교환되어 두 가지 사이에 역인과관계가 존재하는지 판단합니다.

재현성을 높이기 위해, 지역 GWAS 요약-통계 가져오기, SNP 추출, LD 집합, 결과 SNP 추출, 대립유전자 조화, MR 분석, 민감도 분석, 진단 플롯 작성 및 출력 저장에 사용되는 전체 계산 워크플로우가 보충 코드 1로 제공됩니다. 본문은 주요 프로토콜 단계를 설명하며, 보충 파일 1 은 해당 명령어 및 함수 수준 R 구현을 제공합니다.

데이터 출처

요통, 엉덩이 통증, 무릎 통증과 관련된 데이터는 IEU openGWAS 데이터베이스에서 얻으며, 웹사이트는 https://gwas.mrcieu.ac.uk/ 입니다. 이러한 조건을 정의하는 데 사용되는 구체적인 질문은 보충 파일 2에서 확인할 수 있습니다. 요통, 엉덩이 통증, 무릎 통증은 원래 GWAS 데이터셋에서 자기보고한 부위 특이적 통증 항목을 기준으로 정의되었습니다. 요통 데이터셋 GWAS ID 번호는 ebi-a-GCST90018797이며, 표본 크기는 468269이며, 요통 환자 22,413명과 대조군 445,856명, 총 SNP 수는 24174741입니다. 고관절 통증 데이터셋의 GWAS ID 번호는 ebi-a-GCST90013968이며, 표본 크기는 407746이고 총 SNP 수는 11039206입니다. 무릎 통증 데이터셋 GWAS ID 번호는 ukb-b-16254이며, 표본 크기는 461857명이며, 무릎 통증을 가진 98,704명의 환자와 363153명의 대조군, 총 SNP 수는 9851867명입니다. 보행 이상에 대한 GWAS 데이터도 IEU openGWAS 데이터베이스에서 얻었으며, 보행 이상 데이터셋의 GWAS ID 번호는 finn-b-R18 _ABNORMALITI_GAIT_MOBIL이며, 표본 크기는 210717개이며, 1348건의 보행 이상 사례와 209369 대조군을 포함하고 총 1638044개의 SNP 좌위가 확인되었습니다. 자세한 내용은 표 1에 나와 있습니다.

재현 가능한 데이터 수집을 위해 연구자들은 IEU OpenGWAS 데이터베이스에서 위에 나열된 각 GWAS ID를 검색하고, 해당 GWAS 요약-통계 파일을 다운로드한 후 해당 문서를 로컬 텍스트 또는 쉼표 구분된 파일로 저장한 후 R에 가져와야 합니다. 본 연구에서 사용된 계산 워크플로우에서는 소스 파일 형식에 따라 read.table() 또는 read.csv()를 사용하여 로컬 GWAS 요약 통계 파일을 R에 가져왔습니다. 로컬 파일 기반 워크플로우는 동일한 다운로드 GWAS 요약 통계 파일을 고정 열 매핑과 동일한 분석 명령어로 반복적으로 처리할 수 있도록 채택되었습니다. 로컬 파일을 가져오고 처리하는 정확한 R 명령어는 보조 파일 1에 제공되어 있습니다.

SNP의 심사 및 진입

각 노출 표현형과 유의미하게 연관된 SNP는 GWAS 요약 데이터셋에서 P < 5 × 10⁻⁸의 유전체 유의치 값을 사용해 추출하였다. 이 임계값에서 독립적인 SNP가 세 개 미만일 경우, MR 분석을 위한 충분한 유전 기기를 확보하기 위해 P < 5 × 10⁻⁶ 으로 임계값이 완화되었습니다. 연관 불균형은 10,000 kb 창 내에서 r² < 0.001을 사용하여 기기 변수들의 독립성을 보장하였다.

스크립트 수준에서는 노출 GWAS 요약 통계가 먼저 R에 도입되어 미리 정해진 P-값 임계값에 따라 필터링되었습니다. 선택된 노출 SNP는 노출 파일로 저장되어 TwoSampleMR 패키지의 요구사항에 맞게 포맷되었습니다. 필요한 열 매핑에는 SNP 식별자, 효과 추정치, 표준오차, 효과 대립유전자, 기타 대립유전자, 그리고 P 값이 포함되었습니다. 재현 가능한 R 워크플로우에서는 노출 데이터를 read_exposure_data()를 사용하여 가져오고, LD 클럼핑은 클럼프 = TRUE를 설정하여 수행했습니다. 대응하는 R 명령 구조는 다음과 같았습니다: read_exposure_data(filename = "exposure.csv", sep = ",", snp_col = "rsids", beta_col = "beta", se_col = "sebeta", effect_allele_col = "alt", other_allele_col = "ref", pval_col = "pval", clump = TRUE). 열명은 실제로 다운로드된 GWAS 파일에 따라 조정되었습니다.

기기 강도는 F 통계량을 사용하여 평가되었습니다:

F = [(N − k − 1)/k] × [R²/(1 − R²)],

여기서 N은 표본 크기, k는 계측 변수 수, R² 는 SNP로 설명되는 노출 분산의 비율을 나타냅니다. R² 는 다음과 같이 계산되었습니다.

Σ[2 × MAF × (1 − MAF) × β²/(SE² × N)],

여기서 MAF는 소 대립유전자 빈도, β는 대립유전자 효과 추정치, SE는 표준오차입니다.

F 통계량이 10인 SNP> 유지되었습니다. 선천적 해부학적 이상과 체질량지수 등 잠재적 교란 요인과 관련된 SNP를 PhenoScanner V2를 사용하여 확인하여 제거하였습니다.

결과 SNP는 SNP 식별자에 따라 해당 지역 결과 요약-통계 파일과 집합된 노출 SNP를 병합하여 추출했습니다. 결과가 −log10(P)로 보고된 경우, P 값은 P = 10^(−LP)로 변환되었습니다. 추출한 결과 SNP는 read_outcome_data()를 사용하여 임포트되었습니다. 노출 및 결과 데이터셋을 harmonise_data()를 사용하여 효과 대립유전자를 정렬하고 MR 분석에 부적합한 SNP를 제거하였습니다. 대립유전자 방향이 모호한 회문 SNP는 제외하고, mr_keep = TRUE인 SNP는 이후 MR 분석에 포함되었다¹¹. 결과 SNP 추출, P값 변환 및 조화에 관한 상세 스크립트는 보충 파일 1에 제공되어 있습니다.

통계 분석

본 연구에서는 IVW 방법, MR-Egger 회귀법, 가중 중간값 방법, 가중 모드 방법, 단순^{모드 방법 등} 다섯 가지 방법으로 노출 변수와 결과 변수 간의 인과관계를 추정하였다. IVW 방법은 모든 기기 변수가 유효하다고 가정하고 Wald 비율 추정치를 결합하는 원리에 기반하며, 이질성이 있을 때는 무작위 효과 모델을, 이질성이 없으면 고정 효과 모델을 사용하는 MR 분석의 표준 방법으로 간주됩니다. MR-에거 회귀 방법은 잠재적인 다중공선성을 감지하고 회귀¹³에서 인터셉트 항의 존재를 고려합니다. 가중 중앙값 방법은 도구 변수의 50% 이상이 유효한 SNP여야 합니다. 가중 모드 방법은 다른 방법보다 더 작은 샘플 크기를 요구하며, 편향이 적고 타입 I 오류율이 낮아집니다. 단순 모드 방법은 유사한 효과를 가진 SNP를 추정하는 인과 효과가 유사한지 여부에 따라 그룹화할 수 있게 한다¹⁴.

기능 수준에서는 TwoSampleMR 패키지의 mr() 함수를 사용하여 MR 분석이 수행되었습니다. 분석에 사용된 방법 목록은 c("mr_ivw", "mr_egger_regression", "mr_weighted_median", "mr_simple_mode", "mr_weighted_mode"로 지정되었습니다. 오즈비와 95% 신뢰구간은 generate_odds_ratios()를 사용하여 생성되었습니다. 정확한 메서드 목록과 출력 저장 명령어를 포함한 완전한 R 구현은 보충 파일 1에 제공되어 있습니다.

품질 관리 및 민감도 분석

MR 결과의 안정성과 신뢰성을 검증하기 위해 먼저 코크란 Q 검정을 사용하여 SNP의 이질성을 평가하였으며, 코크란 Q 검정이 통계적으로 유의하면 분석된 결과에 유의미한 이질성이 있음을 나타냈다¹⁵. 둘째, MR-Egger 인터셉트 검사를 사용하여 잠재적 수평 다면형성을 평가하였다. 통계적으로 유의미한 MR-Egger 인터셉트는 MR 분석에서 방향성 수평 다면형성의 존재를 나타냅니다. 셋째, 멘델 무작위화 다면형성 잔여합과 이상치(MR-PRESSO)를 적용하여 결과에서 이상치(snp) 존재를 탐색했고, 존재한다면 제거하여 재분석했다¹⁶. 넷째, '리브-원-아웃(leave-one-out)' 방법을 적용하여 결과의 견고성을 검증했으며, 나머지 SNP들의 결합 효과를 하나씩 제외하여 개별 SNP가 노출과 결과 변수 간의 연관성에 미치는 영향을 평가했습니다.¹⁷.

기능 수준에서는 이질성을 mr_heterogeneity(), 방향성 다면형성은 mr_pleiotropy_test(), 단일 SNP 효과는 mr_singlesnp(), mr_leaveoneout()를 사용하여 리브-원아웃 분석이 수행되었습니다. 진단 플롯은 mr_scatter_plot(), mr_funnel_plot(), mr_leaveoneout_plot()를 사용하여 생성되었습니다. MR-PRESSO는 사용 가능한 기기 변수 수가 충분할 때 수행되었으며; SNP 수가 부족해 MR-PRESSO를 수행할 수 없을 경우, 이 과정은 기록되어 보고되었습니다. MR 분석 및 품질 관리 절차는 R 버전 4.3.2와 TwoSampleMR 패키지 버전 0.5.8을 사용하여 수행되었으며, 유의수준은 α = 0.05였습니다. 전체 명령은 보충 파일 1에 제공되어 있습니다.

출력 저장, 품질 검사 및 보고

모든 중간 및 최종 출력은 재현성을 위해 저장되었습니다. 여기에는 선택된 노출 SNP, 집합 노출 도구, 추출된 결과 SNP, 조화된 데이터셋, 보존된 SNP 목록, MR 추정치, 오즈비, 이질성 결과, MR-Egger 인터셉트 결과, MR-PRESSO 결과, 단일 SNP 결과, 생계 결과, 진단 플롯, R 세션 정보가 포함되었습니다. 보고 전에 보유된 SNP는 유전체 전반 유의성, LD 독립성, 대립유전자 조화 상태, 모호한 회문 변이, 기기 강도, 이질성, 방향성 수평 다변성, 이상치 SNP 등을 점검하였습니다. 주요 MR 결과는 95% 신뢰구간과 P값이 포함된 오즈비로 보고되었습니다. 민감도 분석 결과와 진단 도표는 보충 표와 그림에 보고되었습니다. 출력 저장 및 파일 검사 명령은 보조 파일 1에서 제공된다.

결과

요통, 엉덩이 통증, 무릎 통증 등 세 가지 부위특이적 근골격계 통증 표현형과 한 가지 보행 이상 표현형에 대한 GWAS 데이터가 분석되었습니다. 각 노출 인자별 IV에 대한 자세한 정보는 보충 파일 3에서 확인할 수 있습니다.

보행 이상에 대한 부위적 근골격계 통증의 MR 분석

본 연구에서는 연관 불균형과 회문 구조가 있는 SNP를 제외하고 교란변수와 관련된 SNP를 제거했으며, 최종 포함된 SNP는 MR 분석의 도구적 변수로 분석될 것입니다. 모든 SNP의 F-통계량은 10 이상(20.86–91.68)으로, 약한 기기 편향의 가능성이 낮음을 나타냈다. 1차 IVW 분석에서 유전적으로 예측된 요통과 엉덩이 통증은 보행 이상 위험이 높음과 연관되어 있었습니다. 구체적으로, IVW 추정치는 요통(OR = 1.53...

토론

보행 이상은 하체 관절의 부하와 움직임 전략을 변화시키고, 821+_osteoarthritic 과정을 가속화하며, 삶의 질을 크게 저하시키는 기능적 한계를 초래할 수 있습니다. 이러한 임상적 중요성은 보행 기능 장애와 근골격계 통증 간의 상호작용에 대한 관심이 점점 커지고 있습니다. 관찰 연구들은 특정 부위의 근골격계 통증이 보행 장애와 연관될 수 있음을 시사합니다; 예를 들어, 65세 이상 지역사회 거주 노인을 대상으로 한 횡단면 연구에서는 다부위 통증이 보행 속도 느림과 연관되어 있다고 보고했습니다(n = 176)¹⁸. 그러나 무작위 대조 시험(RCT)이 종종 실현 불가능하고, 관찰 결과가 교란 및 역인과에 취약하여 근거 기반 통증 관리 및 재활에 방해가 될 수 있기 때문에 통증-보행 관계의 방향성과 인과관계를 확립하는 것은 여전히 어렵습니다.

최근 여러 연구에서 MR을 활...

재료

이 논문에 사용된 재료 목록

이름	회사	카탈로그 번호	댓글
IEU OpenGWAS 데이터베이스	MRC 통합 역학 부서, 브리스톨 대학교	해당 없음	노출 및 결과 표현형에 대한 GWAS 요약 통계 출처.
MR-PRESSO R 패키지	Marie Verbanck / R 패키지	해당 없음	수평 다면형성의 이상치 검출 및 평가에 사용됩니다.
페노스캐너 V2	케임브리지 대학교	V2	잠재적 교란변수와 관련된 SNP를 식별하는 데 사용되는 웹 자료입니다.
R	R 통계 컴퓨팅 재단	버전 4.3.2	MR 분석 및 품질 관리를 위해 사용되는 통계 컴퓨팅 환경입니다.
TwoSampleMR R 패키지	MRC 통합 역학 부서, 브리스톨 대학교	버전 0.5.8	IVW, MR-Egger, 가중 중앙값, 가중 모드, 단순 모드 MR 분석에 사용됩니다.