Revealing the genetic architectures underlying organ-specific aging based on proteomic data

Domain: neuroscience
Published: 2025-12-11
Editorial outcome: major_revision

Abstract
Organ-specific plasma protein signatures identified via proteomic profiling could be used to quantitatively track organ aging. However, the genetic determinants and molecular mechanisms underlying the organ-specific aging process remain poorly characterized. Here we integrated large-scale plasma proteomic and genomic data from 51,936 UK Biobank participants to uncover the genetic architectures underlying aging across 13 organs. We identified 119 genetic loci associated with organ aging, including 27 shared across multiple organs, and prioritized 554 risk genes involved in organ-relevant biological pathways, such as T cell–mediated immunity in immune aging. Causal inference analyses indicated that accelerated heart and muscle aging increase the risk of heart failure, whereas kidney aging contributes to hypertension. Moreover, smoking initiation was positively linked to the aging of the lung, intestine, kidney, and stomach. These findings establish a genetic foundation for understanding organ-specific aging and provide insights for promoting healthy longevity.

Introduction
Biological aging refers to the gradual and cumulative changes in the biological functions of an organism as time advances1. These changes lead to a decline in physiological capabilities and an increase in disease susceptibility2,3. Investigating the underlying mechanisms of biological aging can enhance our comprehension of the processes contributing to aging-related diseases, thereby informing the development of potential therapeutic interventions.

Most of studies on biological aging focus on whole-body aging indicators, such as healthspan4, lifespan5,6 and longevity7. Given the distinct characteristics of various organs, organ aging, which refers to the biological aging specific to particular organ systems, has increasingly become a research hotspot8–10. Several methods have been established for assessing organ aging. For example, some studies11,12 used MRI-based brain volume and functional connectivity measurements to measure brain aging, but MRI measurements are costly and do not provide molecular insights. Some other studies have used clinical chemistry markers to measure organ aging8,13,14. However, many of these markers have low organ specificity, making them difficult to interpret for organ-specific aging. Moreover, suitable and easily measurable clinical chemistry markers are not yet available for some organs (e.g., intestine and stomach). Therefore, we need to find other indicators for organ aging.

Proteins are the final products of gene expression and the principal functional constituents of cellular and biological processes. Moreover, most drug targets and biomarkers are proteins15,16. In clinical practice, certain plasma proteins have been used to noninvasively evaluate the health status of specific organs, such as alanine transaminase for liver damage17,18. Therefore, finding organ-specific proteins in plasma could track human organ aging. Indeed, recent studies successfully used plasma proteomics data to predict organ aging9,19,20. However, these studies primarily focused on identifying correlations between aging and other phenotypes, without delving into the underlying genetic determinants of organ aging. Moreover, causal associations between organ aging and diseases were not explored.

In this study, using plasma proteomics data from the UK Biobank PharmaProteomics Project (UKB-PPP)21, we trained models to predict organ aging for 13 human organs. Genome-wide association study (GWAS) was then performed to pinpoint the genomic loci associated with organ aging. Mendelian randomization (MR) analysis was further used to explore the causal relationship between aging-related diseases, modifiable lifestyle factors and organ aging. Our results could enhance the comprehension of the genetic basis of organ aging and highlighted the causal relationship between organ aging and other traits.

Results
Overalldesign
TheoveralldesignofthisstudyisshowninFig.1.Briefly,usingplasmaproteomicsdata(2923proteins)from51,936participantsofUKB-PPP21,wemappedputativeorgan-specificplasmaproteomeandusedthemtotrainmodelstopredictorganagingfor13humanorgans(adipose,artery,brain,heart,immune,intestine,kidney,liver,lung,muscle,pancreas,skinandstomach).WethenperformedGWAStoexplorethegeneticmechanismsunderlyingorganaging.Subsequently,weintegratedtranscriptome-wideassociationstudy(TWAS),colocalization,andcausalinferenceapproachestoidentifylikelyeffectorgenes.Furthermore,weconstructedpolygenicriskscorestoassesstheirassociationwithaging-relateddiseases.Additionally,MRanalysiswasemployedtorevealthepotentialcausalitybetweenaging-relateddiseases,modifiablelifestylefactorsandorganaging.

Modelorganagingwithplasmaproteins
Wetrainedagingpredictionmodelsusingproteomicdatameasuredduring2006and2010fromtheUKBiobankPharmaProteomicsProject(UKB-PPP)21for13organs,includingadipose,artery,brain,heart,immunesystem,intestine,kidney,liver,lung,muscle,pancreas,skinandstomach,becauseoftheirrelativelywell-understoodcontributionstoaging-relateddiseases.Theorgan-enrichedproteinsweusedcomprisedtheorgan-enrichedproteinsfromtheHumanProteinAtlas(HPA)project22andtheorgan-specificproteinsidentifiedbyMalmströmetal.‘sstudy23.AsshowninSupplementaryData1,1799proteinswereselectedforthe13organs(rangedfrom34to356).Theseproteinswereenrichedinpathwaysrelatedtoorgan-specificfunctions(SupplementaryData2).Wetrainedagingmodelfor13organsusingorgan-enrichedproteins,anda‘conventional’proteomicagingmodelwastrainedusingallproteins.AsshowninSupplementaryFig.2,amongthethreemachinelearningmodelsused,theelasticnetmodeldemonstratedthebestperformanceonthetestsetsacross11organs.Fortheremainingthreeorgans(kidney,pancreas,andstomach),althoughXGBoostorrandomforestshowedslightlybetterperformanceonthetestdatasetinsomecases,performancediscrepanciesbetweenthetrainingandtestsetswereobservedforthesemodels,indicatingpotentialoverfitting.Incontrast,theelasticnetmodelmaintainedsimilarperformancebetweenthetrainingandtestsets.Therefore,weoptedfortheelasticnetmodelinsubsequentanalyses.TheproteinweightsintheelasticnetmodelareprovidedinSupplementaryData3.Subsequently,wecomparedtheperformanceofthefinaltrainedmodelonthetrainingdatasetandanindependenttestdataset.Themodel’sperformancewassimilaracrossbothdatasets,indicatinggoodgeneralizationability(SupplementaryFig.3).

ThedevelopedmodelswerethenusedtopredictbiologicalagesforallUKB-PPPbaselineparticipants.Similartopreviousstudies24,25,weobservedthatthepredictedagedifference(PAD)—thedifferencebetweenpredictedageandchronologicalage—exhibitednegativecorrelationwithchronologicalage(r rangedfrom−0.90to−0.48,SupplementaryFig.4).Post-predictioncorrectionofthisagebiasproblemiscrucial26.WeadoptthemethodusedbyOhetal.9toaddressthisproblem.ThecorrectedPAD(hereinafterreferredtoas“deltaage”)wasdefinedastheresidualsfromalocallyweightedscatterplotsmoothing(lowess)regressionofpredictedageonchronologicalagemodel.AsshowninSupplementaryFig.4,thismethodsuccessfullycorrectedthedependenceofdeltaageonchronologicalage(rrangedfrom−0.021to0.025).SubsequentanalyseswereperformedafterZ-scorescalingofthedeltaagevalues.

Theassociationbetweenorganaginganddiseasestatus
Wenextexaminedtheassociationsbetweenorganaginganddiseasestatus.DiseasestatuswasdefinedbasedonICD-10diagnoses(UKBiobankField41270).Ninemajordiseasecategorieswereincludedintheranalysis:immunesystemdiseases(ICD-10:D50–D89),circulatorysystemdiseases(I00–I99),musculoskeletalsystemdiseases(M00–M99),nervoussystemdiseases(G00–G99),respiratorysystemdiseases(J00–J99),endocrine,nutritional,andmetabolicdiseases(E00–E90),mentalandbehavioraldisorders(F00–F99),digestivesystemdiseases(K00–K93),anddiseasesoftheskinandsubcutaneoustissue(L00–L99)(SupplementaryData4).Intotal,617specificdiseasesacrossthesecategorieswereincluded.Diseaserecordswerecategorizedintoprevalentdiseases(diagnosedbeforerecruitment)andincidentdiseases(diagnosedafterrecruitment)accordingtothetimingofdiagnosisrelativetorecruitment.

Forprevalentdiseases,weperformedlogisticregressionanalysestoassesstheassociationbetweendiseasestatusanddeltaage.Amongthe126diseasecategory-agingtraitpairs(9diseasecategories×14agingtraits),wedetected125significantpairsofpositivecorrelations.AsshowninFig.2a,somedisease-organcorrelationswereasexpected.Forexample,brainagingexhibitedthestrongestpositivecorrelationwithmentalandbehavioraldisorders.Wealsoobservedsomeotherinterestingcorrelations,forexample,stomachagingwassignificantlycorrelatedwiththestatusofmentalandbehavioraldisorder.Thismightbebecausemultiplecorticalregionscanregulateautonomiccontrolofthestomach,involvingthegut-brainconnection27.Wealsoexploredtherelationshipbetween617diseaseswithinthe9diseasesystemsandorganaging.AfterBonferronicorrection,wedetectedsignificantcorrelationin917disease-organpairs.Foraspecificorgan,thetypeofdiseasethathadthehighestnumberofassociationswithitwasoftentheonethatbelongstothesamesystemastheorgan.Forexample,47.6%(20/42)diseasesassociatedwithacceleratedheartagingwerecirculatorysystemdiseases(Fig.2b).Ontheotherhand,aspecificdiseaseoftenshowedthehighestcorrelationwiththeagingofitsaffectedorgan.Forexample,wefoundthatAlzheimer’sdisease(AD)showedthehighestcorrelationwithbrainaging(Fig.2c).Wealsocategorizedindividualsintoextremeagingtypesusingatwostandarddeviationcutoff(Z-scoreddeltaage>2ase-older,Z-scoreddeltaage<−2ase-younger)foratleastoneagingmodel.Incomparisontoothersubjects,thee-oldercohortexhibitedasignificantlyelevateddiseaseriskwhilethee-youngercohortdemonstratedlowerdiseaserisk(SupplementaryFig.5).

Forincidentdiseases,weappliedCoxproportionalhazardsmodelstoinvestigatetheassociationbetweendeltaageandtheriskoffuturediseaseonset.Amongthe126diseasecategory-agingtraitpairs,wedetected124significantpairsofpositivecorrelations,suggestingthathigherdeltaagewasassociatedwithahigherriskofincidentdiseasestatus.Wealsoexploredtherelationshipbetween617incidentdiseasestatuswithinninesystemsandorganaging.AfterBonferronicorrection,weidentified1526significantcorrelations.Forexample,brainagingwasmostsignificantlyassociatedwithincreasedriskofAlzheimer’sdisease(OR=1.68,p=1.21×10−44)(SupplementaryFig.6).

Validatetheassociationbetweenorganaginganddiseasestatus
Weusedproteomicsdatameasuredafter2019intheUKB-PPP(hereinafterreferredtoasvalidationcohort)studyasanexternalcohorttovalidatetheaging-diseaseassociations.Thisdatasetincluded1463plasmaproteinsfrom1125participants.Individualswhoparticipatedinthevalidationstudywereremovedfromtheoriginalcohort.That’s,thereisnosampleoverlapbetweentheoriginalandvalidationdatasets.

Tomakesurethatthevalidationmodelspreservedmostfeaturesusedintheoriginalmodels,weonlykept8organmodelswithmorethan50%oftheproteinfeaturesavailableinthevalidationdataset(SupplementaryData5).Usingthesefeatures,were-trainedtheorganagingmodelsintheoriginalcohort.AsshowninSupplementaryFig.7,thepredictedagingofthefeature-reducedmodelwashighlycorrelatedwiththatfromthefullfeaturemodel.Inaddition,nosignificantdifferencewasobservedbetweentheweightsfortheproteinsinthefeature-reducedandfullfeaturemodels(SupplementaryData5).Therefore,thefeature-reducedmodelsweresimilartotheoriginalones.Wethenappliedthese8feature-reducedmodelstopredictorganagingforindividualsinthevalidationcohortandanalyzedtheassociationsbetweenpredictedorganaginganddiseaseoutcomes.AsshowninSupplementaryFig.8,wevalidatedtheassociationbetweenorganagingandmultiplediseases.Forexample,weobservedthatadiposeagingwaspositivelyassociatedwithendocrine/metabolicdiseasesintheoriginaldataset(β=0.42,p=1.0×10−111),thisassociationwasalsoobservedinthevalidationdataset(β=0.25,p=1.8×10−6).Thesefindingssupporttheaging–diseaseassociationsweobservedintheoriginalcohort.

GWASoforganaging
Tounderstandthegeneticbasisoforganaging,weperformedGWASonpredicteddeltaage.Thesamplesizeforthe14agingtraits(13organagingand1conventionalaging)rangedfrom19,622to28,151.Aftervariantfilteringandqualitycontrol,upto6,735,132autosomalsingle-nucleotidepolymorphisms(SNPs)wereincludedinouranalysis.

WedividedtheindividualsbasedontheUKBiobankassessmentcentersatwhichparticipantsconsentedintosouthernandnortherncohorts.Weconsideredthesoutherncohort(n=26,293)asthediscoverysamplesandthenortherncohortasthereplicationsamples(n=26,477).Weidentified10,661significantSNPsassociatedwithorganaging(p<5×10−8)inthediscoverysample,98.93%oftheseSNPsalsoexhibitedFDRcorrectedp-value<0.05inthereplicationcohort.Importantly,wefoundthattheeffectsizesoftheseSNPsexhibitedaPearsoncorrelationcoefficientof0.96betweenthediscoverycohortandthereplicationcohort(SupplementaryData6).Finally,wemeta-analysistheGWASresultsfromthetwocohortsandobtained19,767significantSNPsforsubsequentanalysis.Atotalof119loci(definedasa±500-kbregionaroundthegenome-widesignificantSNPs)weredetectedtobeassociatedwithatleastoneagingtrait(Fig.3a,b,SupplementaryFigs.9,10).Wefoundthat92ofthe119loci(77.3%)wereassociatedwithonlyoneorganagingtrait,andlow-to-moderategeneticcorrelationsbetweentheagingof14organswereobserved(averager was0.15)(Fig.3c).

TofurtheridentifypotentialcausalSNPsfororganaging,weutilizedastepwisemodelselectionapproachwithGCTA-COJO28todetectindependentsignals,followedbyFINEMAP29toidentifythemostlikelycausalSNPsforeachsignal.Atotalof1750SNPsachievingeitherconditionalindependenceorahighposteriorprobabilityforcausality(log BF>2)wereretainedandreferredas‘fine-mappedSNPs’forsubsequentanalysis.Weidentified7loci-thathadnomorethan5variantswithinthe95%credibleset(Fig.3d),whichcanbeconsideredashighlycrediblecausalSNPs.

Only78fine-mappedSNPs(5.21%)werelocatedwithinexonofprotein-codinggenes,including47nonsynonymousSNPs.Tofurtherinvestigatethefunctionoffine-mappedSNPs,weannotatedtheseSNPsusingpublicdatabases.AsshowninSupplementaryData7,7SNPswereannotatedashavingclinicalandpathologicalsignificanceintheClinVar30database,165SNPswerepredictedtobedeleteriousto-genefunctionbyCADD31,DANN32,orFATHMM-MKL32,1174SNPswerepredictedaslikelytoaffectregulomebindingbyRegulomeDB33,9SNPswerepredictedtoaffectgenesplicingbyregSNPs-intron34,and30SNPswerepredictedasdamagingexonicvariantswiththepotentialofaffectingproteinfunctionbySIFT/PROVEAN35,36.

Queryingfine-mappedSNPsfromexistingGWAS
Toinvestigatewhetherthefine-mappedSNPsassociatedtoorganagingexhibitknownassociationswithotherphenotypes,wequeriedtheseSNPsintheGWAS-catalogdatabase37.Theresultsshowedthat15.8%(276/1750)ofthefine-mappedSNPshavebeenreportedtobeassociatedwithothertraits(SupplementaryData8).Wefoundthatfine-mappedSNPswereannotatednotonlytothecorrespondingorganphenotypebutalsotophenotypesofotherorgans.Forexample,Thegenomiclocusexhibitingthemostsignificantassociationsignalwithlungaging(chr10:80716529–82718385)hasbeenreportedtobeassociatedwithbothlung-relatedtraits,suchassuchastheforcedexpiratoryvolumeinonesecondtoforcedvitalcapacityratio(FEV1/FVC)38andCOVID-1939,butalsowithotherorganrelatedtraits,suchasrheumatoidarthritis40andheight41.Alocusassociatedwithbrainaging(chr1:149714166–151444734)hasbeenreportedtobe-linkedtoschizophrenia42andphenotypesofotherorgans,includingbloodproteinlevels43,plateletcount44.

Geneprioritizationfororganaging
TolinkGWAS SNPstorelevantprotein-codinggenes,weemployedgene-levelandSNP-levelprioritizationmethods.Thegene-levelanalysesincludedTWAS,colocalization(usingtheRpackagecoloc45),andsummary-data-basedMendelianRandomization(SMR)46.IntheSNP-levelanalyses,genesharboringfine-mappedSNPs-withfunctionalannotationsupportwereincludedintheprioritizedgenes.Finally,weidentified554protein-codinggenesassociatedwithorganaging(Fig.4a,SupplementaryData9).

Theimmunesystemharboredthegreatestnumberofriskgenes(n=91),whilemusclehadonly10riskgenes(Fig.4b).Amongthe554riskgenes,only37wereassociatedwithmultipleagingtraits.FUT2,whichwasariskgeneofintestine,pancreas,andstomachaging,wasreportedtoberelatedtothefunctionsoftheseorgans47–49.AnothergenePLG,whichwasreportedtobeassociatedwithcardiovascularfunction50andliverregeneration51,wasidentifiedasariskgeneforarteryandliveraging.

ABObloodtypeisassociatedwithaginginmultipleorgans
Amongthe37riskgenesassociatedwithmultipleagingtraits,thegeneABO,whichdeterminesbloodtype,wasidentifiedasariskgeneforadipose,artery,immune,intestine,pancreas,andlungaging.IntheUKBiobank,bloodtypewasinferredbasedontwoSNPs:rs8176719andrs8176746.IntheGWASresults(Fig.4c),Talleleofrs8176719waspositivelyassociatedwiththeagingofartery(β=0.05,p=1.4×10−8),immune(β=0.04,p=7.7×10−10),intestine(β=0.06,p=9.6×10−12)butnegativelyassociatedwithpancreasaging(β=−0.07,p=2.7×10−15).Galleleofrs8176746waspositivelyassociatedwithpancreasaging(β=0.17,p=1.5×10−23).Asforlungaging,nosignificantresultswerefoundforbothSNPs,butsuggestivesignificantsignalwasdetectedforrs8176719(βforTallele=−0.03,p=2.1×10−6).

Wefurthercomparedthedeltaageacrossdifferentbloodtypes.WefoundthatindividualswithbloodOtypeexhibitedthehighestarterialandintestinalage.ThosewithbloodtypeBdisplayedthelowestlungandpancreasage.Meanwhile,individualswithbloodtypeAhavethelowestlungagecoupledwiththehighestpancreasage(SupplementaryFig.11).Thesefindingssuggestedthatpeoplewithcorrespondingbloodtypesshouldpayspecialattentiontothehealthstatusofrelatedorgans.

Functionalanalysisoftheorganaginggenes
Pathwayenrichmentofriskgenes.Toannotaterelevantbiologicalpathwaysunderlyingeachorgansystem,weperformedpathwayenrichmentanalysisfortheorganaginggenes.Wefoundthatsomeorgan-enrichedpathwaysshowedstrongorganspecificity,reflectingtheuniquefunctionsoftheseorgans.Forexample,thetop-rankedenrichedpathwaysforadiposeaginggenesincludedfattyacidbeta-oxidation;themostsignificantlyenrichedpathwayforheartaginggenesincludedsarcomereorganizationandmyofibrilassembly;themostsignificantlyenrichedpathwaysforimmuneaginggenesincludedTcellmediatedimmunity.Ontheotherhand,someorgan-enrichedpathwaysexhibitbroaderfunctions.Forinstance,themostsignificantlyenrichedpathwayforintestineaginggeneswasepidermalgrowthfactorstimulus,whichwasrelatedtoprocessessuchascelldifferentiationandregulationofinflammation(Fig.5a).

Annotationofriskgenesinmousephenotypes.TheMouseGenomeInformatics(MGI)52databaseprovidesphenotypicinformationresultingfromgenedeletionsormutationsinmice.Thesephenotypeswereclassifiedintodifferentcategoriesaccordingtotheirrelatedorganorbiologicalprocess.Toinvestigatethefunctionalrelevanceoforganaginggenes,wetestedwhetherriskgeneswereenrichedinspecificphenotypiccategories.AsshowninFig.5b,wefoundthatadiposeaginggenesweresignificantlyenrichedinadiposetissuephenotypes,suggestingthatdisruptionofthesegenesinmiceislikelytocauseabnormalitiesinadiposetissue.Heartaginggeneswereenrichedinmuscleandcardiovascularsystemphenotypes,whileimmuneaginggenesshowedstrongenrichmentinimmunesystemandhematopoietic-relatedphenotypes.

Gene-drug-diseasesnetwork.Toinvestigatethepotentialoforganaginggenesasdrugtargets,weusedGREP53toperformdrugtargetenrichmentanalysis.GREPreportedannotateddrugsofinputgenesandthecorrespondingenricheddiseasecategories,therebyconstructingagene–drug–diseasenetworkofpotentiallyrepositionabledrugs.AsshowninFig.5c,twogenesrelatedtobrainaging(APH1A,SV2A)weretargetsof2drugs,whichwereusedtotreatextrapyramidalandmovementdisorders(ICD-10codeG20-G26).PLGisariskgeneassociatedwithartery,liver,andconventionalaging,andwasusedtodevelopdrugsfortreatingdiseasesoftheear,gallbladder,andpuerperium,indicatingpotentialfordrugrepurposing.

Organagingriskprediction
Toinvestigatewhetherthegeneticorganagingcouldpredictdiseasestatus,weperformedpolygenicriskscore(PRS)analysis.PRSwascalculatedinUKBparticipantsthatarenotpresentinUKB-PPPtoensurethatthereisnooverlapinparticipantsbetweentheorganagingsummarystatisticsandthetestedindividuals.Wedefinedthetop5%ofeachorgan’sPRSashighagingriskindividualsandthebottom5%aslowagingriskindividuals.Subsequently,weusedlogisticregressiontocalculatetheassociationofhigh-riskandlow-riskindividualsandninediseasecategories.Inhighagingriskindividuals,mostdisease-organpairsshowedapositivecorrelation(98/126),with9remainingsignificantlypositiveaftermultipletestingcorrection.Intheanalysisoflow-riskindividuals,mostdisease-organpairsshowedanegativecorrelation(90/126),with4remainingsignificantlynegativeaftermultipletestingcorrection(Fig.6a).Wealsoexploredtherelationshipbetween617diseaseswithin9diseasecategoriesandagingrisk.AfterBonferronicorrection,62significantdisease-organpairsweredetected,and85.48%(53/62)ofthemwerepositiveassociations(SupplementaryData10).

Mendelianrandomization(MR)
Toinvestigatethecausaleffectsoforganagingondiseasestatusandexplorepotentiallifestyleinterventionsfororganaging,weperformedtwo-sampleMR54analysis.TheassociationoftheIVsandorganagingwasderivedfromtheGWASresultsinthisstudy.

Causaleffectsoforganagingondiseasestatus.Weutilizedtwo-sampleMR54toexaminethecausaleffectsoforganagingon18commondiseasesinvolvingmultipleorgans.TheoutcomedatawasobtainedfromtheFinnGen55database(detailinSupplementaryData11).Weidentified31potentialcausaleffectsoforganagingondiseases(p<0.05),withthreepairsremainingsignificantafterBonferronicorrection(Fig.6b).Leave-one-outanalysesshowedthatnosingleSNPdrovethecausalestimates(SupplementaryFig.12).NosignificantpleiotropywasdetectedusingtheMR-PRESSOglobaltestorMR-Eggerintercepttest.AsshowninFig.6b,bothheartandmuscleagingareassociatedwithincreasedriskofheartfailure(p=1.84×10−5and2.14×10−5,respectively);kidneyagingispositivelyassociatedwiththeriskofhypertension(OR(95%CI)=1.06(1.03,1.09),p=1.71×10−4).

Tovalidatetheidentifiedsignificantcausaleffectoforganagingonheartfailureandhypertension,were-performedMRanalysisusingoutcomeGWASsummarystatisticsfromtwoadditionalcohorts.TheobesitydatasetwasreleasedbyShah et al.56 and the hypertension dataset was published by Singh et al.57. The exposure data was the same, that’s, the IVs we used for validation the causal effects were the same specific SNPs in the original MR analysis. The results (Supplementary Fig. 13) are consistent with those from the FinnGen datasets, confirming that both heart and muscle aging were linked to an elevated risk of heart failure (p=7.28×10–7 and 9.62×10–5, respectively), while kidney aging correlated with an increased risk of hypertension (p=1.32×10–4). Leave-one-out analyses showed that no single SNP drove the causal estimates (Supplementary Fig. 14).

Potentiallifestyleinterventionsfororganaging.WealsousedMRtoexplorepotentiallifestyleinterventionsfororganagingusinglifestylesasexposureandorganagingasoutcome.Thelifestylesincludeageofinitiationofregularsmoking(AgeSmk),whetheranindividualhadeversmokedregularly(SmkInit),cigarettesperday(CigDay),smokingcessation(SmkCes)anddrinkspersweek(DrnkWk)58.AsshowninFigs.6c,5significantassociationswereidentifiedaftermultipletestingcorrection.Forexample,SmkInitispositivelyassociatedwithlungaging(β(95%CI)=0.27(0.20,0.34),p=2.38×10−8),intestineaging(β(95%CI)=0.23(0.13,0.33),p=1.98×10−6),kidneyaging(β(95%CI)=0.11(0.05,0.16),p=9.65×10−5)andstomachaging(β(95%CI)=0.26(0.17,0.35),p=2.41×10−7).Theeffectofthesmokingphenotypeonacceleratinglungagingwasanexpectedoutcome.DrnkWkexhibitedasuggestivepositiveassociationwithbrainaging(β(95%CI)=0.19(0.03,0.34),p=0.018)andconventionalaging(β(95%CI)=0.21(0.06,0.36),p=5.11×10−3)(Fig.6c).

Methods
Theoriginalandvalidationcohort
TheoriginalcohortoftheproteomicsdatawasobtainedfromtheinitialrecruitmentphaseoftheUKBiobank(measuredin2006–2010,including2923proteinsmeasuredfrom53,021subjects).Anotherdatasetfromthefollow-upstudyofUKBiobank(measuredin2019andlater,including1463proteinsmeasuredfrom1125subjects)wasusedtovalidatetheaging-diseaseassociationsweobservedfromtheoriginalcohort.Toavoidsampleoverlap,weremoved1085individualsfromtheoriginalcohortwhohadparticipatedinthevalidationstudy.Therefore,thesamplesizeforouroriginalstudywas51,936.

Organ-enrichedproteins
Theorgan-enrichedproteinsweusedwerethecombinedsetoforgan-enrichedproteinsfromtheHumanProteinAtlas(HPA)project22andtheorgan-specificproteinsidentifiedbyMalmströmetal.’sstudy23.Briefly,intheHPAproject,thetissuespecificityofagenewasbasedonmRNAexpressionlevelsintheconsensusdataset,includingtheGTEx67datasetsandHPAdatasetsfrom1206celllines,40humantissuesand18immunecelltypes.Onlyproteinswithexpressionlevelsinaparticulartissue/region/celltypeatleastfourtimesanyothertissue/region/celltypewereconsideredasorgan-enrichedproteins.Asfororgan-specificproteinsfromMalmströmetal.’sstudy23,theydevelopedamass-spectrometry-basedproteomicsstrategytoinferthetissueoriginofproteinsdetectedinhumanplasma.Onlyproteinsidentifiedasoriginatingsolelyfromaspecificorganwereconsideredasorgan-enrichedproteinsforthatorgan.IntheHPAdataset,proteinsoriginatingsolelyfromthelymphoidorbonemarrowwereclassifiedasimmunesystem-enrichedproteins.InMalmströmetal.’sstudy23,proteinsoriginatingsolelyfromthespleenorbonemarrowwereclassifiedasimmunesystem-enrichedproteins.Detailsoftheorgan-enrichedproteinswereprovidedinSupplementaryData1.

Imputationofplasmaproteomicdata
Intheoriginalstudy,foreachorgan,weextractedproteomicsdataoforgan-specificgenesandexcludedindividualswithamissingrategreaterthan10%.Thismeansthatfororganswithfewerthan10specificgenes,weselectedindividualswithoutanymissingvalues.Next,weusedtheKNNalgorithmtoimputetheproteomicsdataforeachorgan.ImputationwasperformedusingRpackageimpute(v3.3.6)68.

Agingmodeltraining
Wetrainedmodelsusingproteomicdatafromhealthyindividuals(definedashavingnoICD-10diseaserecords,n=7217).Foreachmodel,onlyorgan-enrichedproteinswereincludedasinputfeatures,andsexwasincludedasacovariate.Weemployednestedcross-validationtoassessmodelperformanceandoptimizehyperparameters.Theentiredatasetwasdividedintofiveouterfolds,withonefoldheldoutasthetestsetineachiteration.Withinthetrainingdataofeachouterfold,weusedafive-foldinnercross-validationtoperformhyperparametertuningusingGridSearchCVfunctioninscikit-learn(v1.4.1)69.Foreachsetofhyperparameters,themodelwastrainedonfourinnerfoldsandevaluatedontheremainingfold,andthebest-performinghyperparameterswereselected.Thecompletenestedcross-validationprocedureisshowninSupplementaryFig.1.Weselectedthreemachinelearningalgorithms:elasticnet,XGBoost,andrandomforest.Thefollowingparametersweretunedformodeltraining:

Elasticnet.Thealphaparameter,whichcontrolstheoverallstrengthofregularization,wastestedfrom0to1inincrementsof0.01.Similarly,thel1_ratio,whichdeterminesthebalancebetweenL1andL2penalties,wasalsovariedfrom0to1instepsof0.01.

XGBoost.Thelearningratewastestedatvaluesof0.01,0.1,and0.3,whilethenumberofestimatorsrangedfrom100to500instepsof100.

Randomforest.Thenumberofestimatorsrangedfrom100to500instepsof100.

Modeldevelopmentandevaluationswereimplementedthroughxgboost(v2.1.1)andscikit-learn(v1.4.1)underPython(v3.10.10).Weusedthedatafrom1085individualswhoparticipatethefollow-upstudyoftheUKBB(measured2019+)astheindependenttestdataset.Toavoidsampleoverlap,individualsfromthetestdatasetwereexcludedfromtheoriginalcohortattheverybeginning.Weremovedindividualswithdiagnosesfromtheninemajordiseasecategories(detailsinSupplementaryData4),resultinginanindependenttestdatasetof622participants.Thesexandagedistributionsofthetestdatasetwerecomparabletothoseobservedinthetrainingdataset(SupplementaryFig.16).

Deltaagecalculation
Previousstudies24,25haveidentifiedagebiasasacommonprobleminagingpredictionmodels,wherethepredictedagedifference(PAD)—thedifferencebetweenpredictedageandchronologicalage—tendstobenegativelycorrelatedwithchronologicalage.Post-predictioncorrectionofthisbiasiscrucial26.WeadoptthemethodusedbyOhetal.9toaddressthisproblem.WecalculatedthecorrectedPAD(deltaage)astheresidualsfromalocallyweightedscatterplotsmoothing(lowess)regressionofpredictedageonchronologicalagemodel.Hereweusedlowessregressioninplaceofasimplelinearregressionbecauseextensiveevidenceshowedthattheplasmaproteomechangesnon-linearlywithage70.Lowessregressionwasimplementedusingstatsmodels(v0.14.5)underPython(v3.10.10).

Estimationofthecorrelationbetweendeltaageanddiseasestatusorchronologicalage
Diseaserecordswereclassifiedintotwocategories:prevalentdiseases(diagnosedpriortorecruitment)andincidentdiseases(diagnosedafterrecruitment).Logisticregressionwasusedtoestimatethecorrelationofdeltaageandprevalentdiseasestatus,withage,sexandPC1-5ascovariates.Coxproportionalhazardsmodelswereusedtoassesstheassociationbetweenorganagegapsandincidentdiseasestatus,modelsweretestedusingthelifelinespythonpackage71withage,sexandPC1-5ascovariates.Individualswiththediseasewerecodedas1andthosewithoutthediseaseas0.Weusedlinearregressiontoestimatethecorrelationofdeltaageandchronologicalage,withsexandPC1-5ascovariates.Logisticregressionwasimplementedusingstatsmodels(v0.14.5),Coxproportionalhazardsregressionwasimplementedusinglifelines(v0.30.0),andthePythonversionusedwasv3.10.10.

Validationofaging-diseasesassociationinotherdatasets
Wechosetheproteomicsdatafromthefollow-upstudyoftheUKBB(measured2019+)asthevalidationdataset.Thisdatasetincluded1463plasmaproteinsfrom1125participants.Tomakesurethatthevalidationmodelspreservedmostfeaturesusedintheoriginalmodels,weonlykeptorganmodelswithmorethan50%oftheproteinfeaturesavailableinthevalidationdataset(SupplementaryData5).Usingthesefeatures,were-trainedtheorganagingmodelsintheoriginalcohort.Thesefeature-reducedmodelswerethenusedfororganagingpredictioninthevalidationdatasets,andaging-diseaseassociationwerealsocalculated.

GWAS
WeperformedGWASof14agingtraitsafterremovingindividualswithZ-scoreddeltaageabovemean+3*SDorbelowmean−3*SD.Analyseswererestrictedtovariantswithminorallelefrequency>0.01,genotypehardcallrate>0.98,HardyWeinbergp>1×10−7,andpresentintheHaplotypeReferenceConsortium(HRC)panel72.Participantswereexcludedfromanalysisiftheymetanyofthefollowingcriteria:(1)mismatchbetweenself-reportedsexandgeneticsex,(2)genotypingcallrate<0.95,(3)wereoutliersforheterozygosity(deviatemorethan3standarddeviationsfromtheheterozygosityratemean),or(4)relatedindividuals(only1individualinpairwithapihatgreaterthan0.044wasremained).Finally,6,735,132autosomalvariantsand33,243samplesinUKB-PPPwereremained.GWASwasperformedusingamixedlinearmodelimplementedinBOLT-LMM(v2.4.1)73,adjustingforage,sex,age-squared,age-sexinteraction,age-squared-sexinteraction,thetop40geneticPCs,bloodpressure,andBMI.WedividedtheUKBBdatasetintotwogroups,southernandnorthern,basedontheenrollmentlocationoftheUKBiobankassessmentcenter(FieldID54).TheseassessmentcentersaredistributedacrossEngland,Scotland,andWales.AccordingtothestandardsoftheUKOfficeforNationalStatistics(ONS)74,EnglandwasfurtherdividedintoSouthernEngland,NorthernEngland,andtheMidlands.Basedonthegeographicallocationsoftheseregions,individualsfromSouthernEngland,theMidlands,andWaleswereclassifiedastheSoutherncohort(n=26,293),whilethosefromNorthernEnglandandScotlandwerecategorizedastheNortherncohort(n=26,477).Meta-analysesforGWASresultsobtainedfromdiscoveryandreplicationsampleswasperformedusingMETAL75.TheSNPheritabilityandgeneticcorrelationwerecalculatedusingLDSC(v1.0.1)76. The genomic coordinate is according to Genome Reference Consortium human build 37 (GRCh37) in all analyses.

Fine-mappingSNPs
First,SNPsinsignificantwere-definedasbeingconditionallyindependentusingGCTA-COJO(v1.95.0)28. Subsequently, we utilized FINEMAP(v1.4.2)29 software to calculate posterior probabilities for all SNPs within 500 kb of each conditionally independent SNP.

Functionalannotationofthefine-mappedSNPs
Weannotatedthefunctionoffine-mappedSNPsusingClinVar30,CADD31,DANN32,FATHMM-MKL77,RegulomeDB33,regSNPs-intron34,andSIFT/PROVEAN35,36. SNPs supported by at least one database annotation were considered as functional SNPs. Detail annotation criteria are as follows: 1) SNPs annotated as ‘pathogenic’ or ‘likely pathogenic’ in ClinVar30 were considered to be having clinical significance; 2) SNPs with score ≥12.37 in CADD31, or DANN32 score ≥0.9, or FATHMM-MKL77 score ≥0.7 were considered as deleterious to gene function; 3) SNPs with RegulomeDB33 rank ≤2 were considered as affecting regulome binding; 4) SNPs annotated as “damaging” in regSNPs-intron34 were considered as affecting splicing; 5) SNPs with SIFT39 score <0.05 or PROVEAN40 score <−2.5 were considered as damaging exonic variants.

Geneprioritization
TheeQTLdataforTWAS,COLOC,andSMRanalyseswereobtainedfromGTExV867. For the 13 types of organs, only the data from the corresponding tissues were utilized, and for conventional aging, all tissues were used. Then, the best results across tissues were selected. For TWAS and SMR, we retained genes that remained significant after Bonferroni correction. For COLOC, we retained genes with posterior probability of hypothesis 4 (PPH4) greater than 0.9. TWAS was performed based on imputed cis-regulated gene expression using FUSION with default settings46,78,79. Pre-computed gene expression weights of 49 tissues from GTExv8 were used as downloaded from the FUSION website (http://gusevlab.org/projects/fusion/). The Bonferroni correction threshold was defined as 0.05 divided by the total number of scored genes in that organ. COLOC was performed using R package coloc(v5.1.0.1)45 with default settings. SMR was performed using SMR(v1.4.0)46 with default settings. Genes containing fine-mapped SNPs with functional annotation support (as described in the functional annotation of the fine-mapped SNPs part) were also included in the prioritized genes. SNPs were assigned to their corresponding genes using ANNOVAR80.

Comparisonofdeltaageacrossdifferentbloodtypes
Wefirstlycalculatedtheresidualsofdeltaageafterregressingitonage,sex,age-squared,age-sexinteraction,age-squared-sexinteraction,thetop40geneticprincipalcomponents(PCs),bloodpressure,andBMI.Theresidualsrepresenttheunexplainedvariationindeltaageafteraccountingforthesepotentialconfounders.Wethencomparedtheresidualdeltaagevaluesacrossdifferentbloodtypegroups.Thestatisticalsignificanceofbetween-groupdifferenceswasestimatedusingatwo-sidedrank-sumtest.

Pathwayenrichment
PathwayenrichmentwasperformedusingRpackageclusterProfiler(v4.16.0)81,wheregenesymbolswereconvertedtoEntrezIDsbeforeconductingtheenrichmentanalysis.

EnrichmentofprioritizedgenesinmousephenotypecategoriesfromMGI
MGI52isaninternationaldatabaseresourceforthelaboratorymouse,providingphenotypicinformationderivedfromgenedeletionsormutationsinmice.Thesephenotypesarecategorizedbytheirassociatedorganorbiologicalprocess.WeperformedFisher’sexacttesttotestwhethertheorganaginggeneswereenrichedinspecificcategories.

Gene-drug-diseasesnetwork
DrugtargetenrichmentanalysiswasperformedusingGREP(v1.0.0)53. Foreachorgan,theinputwasriskgenelist,GREPreportedthedrugsannotatedtotheseriskgenesandthetypesofdiseasestheyareenrichedin.

PRS
WeusedPRS-CS(v1.1.0)82tobuildagingriskpredictionmodels.PRS-CSisaPython-basedcommand-line-toolthatusesGWASsummarystatisticsandanexternalLDreferencepaneltoinfertheposterioreffectsofSNPimpactlevelsundercontinuousshrinkage(CS)priorconditions.Inthemodel,weincluded1,161,863SNPsfromHapmap3. For the reference LD panel, we used data constructed from the UK Biobank (https://www.dropbox.com/s/t9opx2ty6ucrpib/ldblk_ukbb_eur.tar.gz?dl=0). We randomly selected 50,000 unrelated samples from the UKBB dataset as the target dataset. For other parameters, we opted for the default values. The output of PRS-CS contained the chromosome, rs ID, base position, A1, A2, and posterior effect size estimate for each SNP. Individual-level polygenic scores were produced by concatenating output files from all chromosomes and then using PLINK’s (v1.9.0) --score command (https://www.cog-genomics.org/plink/1.9/score)83.

MR
We collected GWAS summary statistics of 18 diseases from the Finngen55 database, and GWAS summary statistics of 5 smoking and drinking phenotypes from a study based on non-UKBB populations84. The association of the IVs and organ aging was derived from the GWAS results in this study.

InstrumentVariant(IV)selection.WeusedtheclumpingalgorithminPLINK(v1.9.0)(https://www.cog-genomics.org/plink/1.9/)83 to select independent SNPs for each exposure (r2 threshold=0.001, window size=1Mbandp<1×10−5).The1000GEuropeandata(phase3)wereusedasthereferenceforLDestimation85. Secondly, if one selected SNP in the first step was not present in outcome data, we used its proxy SNP with r2>0.8 instead. The proxy SNP must also be associated with the exposure (p<1×10−5). If more than one proxy SNPs were available, we chose the one with the maximum r2 and minimum p value associated with the exposure.

IVqualitycontrol.FortheIVs,threekeyassumptionsmusthold:1)theselectedIVsmustbeassociatedwiththeexposure(relevanceassumption);2)theselectedIVsarenotassociatedwithpotentialconfounders(independenceassumption)and3)theIVsaffecttheoutcomeonlythroughtheireffectontheexposure(exclusionrestrictionassumption).WeusedRadialMRpackage(v0.1.0)86underR(v4.3.1)toremovepleiotropicSNPs.TheremainingSNPswereusedtoperformMRanalysis.WeusedMRSteigerfilteringtocheckwhethertheMRanalysisestimatesassessedthetruecausaldirection87.

Bidirectionaltwo-sampleMRanalyses.WeselectedIVWregressionwithmultiplicativerandomeffectsmodelastheprimarymethodforcausalinference88. That’s, the evidence of statistically significant causal effect was based on the IVW p-value. To complement and enhance the reliability of the results, we conducted four other MR methods. MR-Egger method estimates the causal effect through the slope coefficient of the Egger regression, which provides a more robust estimate even if none of the IVs are invalid89. Weighted median method can even protect against up to 50% of invalid IVs90. Weighted mode method provides consistent estimates when the relaxed IV assumption has less bias and a lower type-I error rate91. MR-RAPS accounts for systematic and idiosyncratic pleiotropy and can provide a robust inference for MR analysis with many weak instruments92. All these methods were implemented in the TwoSampleMR93 package (v0.6.8) under R (v4.3.1).

Sensitivity analysis. For the significant MR results after Bonferroni correction, we further performed sensitivity analysis. First, we performed leave-one-out analysis to check whether the causal association was obviously driven by a single SNP. Second, we conducted MR-PRESSO to detect the presence of horizontal pleiotropy (p<0.05)94. Third, we executed MR-Egger regression to examine the potential bias of directional pleiotropy89. The intercept in the Egger regression indicates the mean pleiotropic effect of all genetic variants, which is interpreted as evidence of directional pleiotropy when the value differs from zero (p<0.05). MR-PRESSO was implemented in the MRPRESSO package (v1.0) under R (v4.3.1). MR-Egger was implemented in the TwoSampleMR93 package (v0.6.8) under R (v4.3.1).

DatasetsusedforthevalidationofMRresults
To validate the identified significant causal effect of organ aging on heart failure and hypertension, we obtained outcome GWAS summary data from two other datasets to re-perform MR analysis. The heart failure data including 583,167 non-UKBB European individuals was from the study by Shah et al.56, and the hypertension data including 10,772 African individuals was from the study by Singh et al.57. The summary data for exposure were the same as those used in the primary analysis. That’s, The IVs we used for validation the causal effects were the same specific SNPs in the original MR analysis.

Statistics&reproducibility
Nostatisticalmethodwasusedtopredeterminesamplesize.Duringqualitycontrolofproteomicandgenomicdata,individualsnotmeetingpredefinedqualitycriteriawereexcludedfromtheanalyses.Thenestedcross-validationprocedureusedinmodeltraininginvolvedrandomsplittingofdata;therandomnesswascontrolledbysettingafixedrandomseedtoensurereproducibility.

Ethics
TheuseofUKBiobank(UKBB)datawasapprovedbytheUKBBaccordingtoitsestablishedaccessprocedures.TheUKBBhasethicalapprovalfromtheNorthWestMulti-centreResearchEthicsCommitteeasaResearchTissueBank.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "In this study, the authors constructed models to predict organ-specific aging using proteomics data, and conducted a thorough investigation searching for genetic determinants of organ-specific aging in large populations. Extensive functional annotations were performed to explore effector genes of the identified loci. They also construct polygenic risk scores to assess associations with aging-related diseases for clinical implications. This study is well designed and provides very interesting results of genetic architectures for organ aging, which provide new insights into the understanding of the aging process and facilitate the discovery of anti-aging interventions. The methodology are sound to support the conclusions. The manuscript is well-written. Overall, this is a very nice piece of work, which deserves to be published. I have only a few minor points to potentially improve this work.\nComments:\n1. The authors looked up whether the significant loci associated with corresponding organ-associated phenotypes in previous GWAS studies. It would be better to check the previous association results for phenotypes related to other organs, since \"no organ system is an island\".\n2. When comparing with previous GWASs, the authors should present the results of overlapping findings.\n3. Methods: The authors should provide the weights of proteins in each model to help readers understand the importance or contribution of each protein in the model’s predictions.\n4. Specify the method used to assess the relationship between delta ages and diseases, which should be explicitly described in the text.\n5. The authors should make the GWAS summary data publicly available.\n6. In line 277, \"supplementary Table 17\" should be corrected to \"supplementary Table 16\".\n7. The X-axis title is missing in Figure 5A.", "reviewer_2": "Zhu and colleagues used Olink plasma proteomics data from the UK Biobank to study the genetic architecture of the multi-organ proteomics age derived from machine learning models. Proteomics age has been an increasingly used biomarker in aging research. The authors used three different machine-learning models to derive the organ-specific age. Then, they used the Z-scores aging biomarkers as phenotypes to run a GWAS and several post-GWAS analyses. The authors made a great effort to analyze the data extensively. The topic is needed in a timely manner and may provide additional insights on previous research endeavors. However, the rigor and interpretation of these results need to be further scrutinized. Below are some of my concerns and comments.\n1) Writing: Overall, the writing and flow of the current paper need improvement. Especially for the Abstract, the main messages are obscure, and it seems more like a Method section, not a highlight of the critical findings in this study. For example, the authors state: \"Furthermore, we construct polygenic risk scores to assess their association with aging-related diseases. Subsequent Mendelian randomization analysis revealed the potential causal effects of organ aging on disease risk status\". I don't think these are the main messages. Instead, they are the specific methods and analyses the authors used in this study. Similarly, for the rest of the manuscript, a clear flow should guide the organization of the writing, including the intuition for the genetic analyses you chose to run.\n2) ML model for age prediction: Regarding the age prediction model, the authors should carefully scrutinize the model's overfitting/generalizability with the ML models you employed on the proteomics data. In the previous papers (https://www.nature.com/articles/s41586-023-06802-1; https://www.nature.com/articles/s41591-024-03164-7), the authors found evident model overfitting and used feature selection or regularization terms to battle this phenomenon. In the current writing, the author only showed the Pearson's r between the predicted age and chronological age in Supplementary eFigure S2 (not sure whether this is for test data and the training or the entire population). The authors should report Pearson's r and other loss metrics for the training/validation/test results to show no severe overfitting. In S Fig. 3, the authors associated the chronological age with the delta age and showed no significant associations. This is very surprising, as age bias (https://www.sciencedirect.com/science/article/pii/S2213158223000086) is very common in any age prediction model. In addition, why did the author use logistic regression to associate two quantitative variables?\n3) Genetics: For GWAS, the QC'ed number of valid SNPs (13206758) is much larger than that of the previous GWAS, using imputed data from the UK Biobank. Can authors state the discrepancy?\nSee some references from previous UKBB GWAS:\nhttps://www.nature.com/articles/s41586-018-0571-7#Sec8\nhttps://www.nature.com/articles/s41467-020-17368-1?fromPaywallRec=false#Sec9\nhttps://www.nature.com/articles/s41588-019-0516-6\nThe covariates included in the GWAS are insufficient, especially since proteomics aging has been shown to be non-linear. I suggest adjusting for the effects of age, sex, age-squared, age-sex interaction, age-squared-sex interaction, and the top 40 genetic PCs, as well as blood pressure, BMI, etc. In addition, the authors should run the GWAS using a mixed linear model instead to better combat aging population stratification. The Manhattan and QQ plots are missing, which is critical to check the quality of the GWAS.", "reviewer_3": "This submitted manuscript is an interesting work that mainly explores the interplay between aging and disease, via computational analyses on genetic variants, genes and proteins in a tissue-specific manner as well as causal inference. Nonetheless, I have the following major concerns before any further evaluation of the scientific significance as well as the fidelity of methods and results.\n(1) Lack of validation or replication procedures in the study design.\nAll results, including the proteomics imputation for aging rates, then GWAS and PRS, as well as the MR analysis and thereafter, were all discoveries without replication. The authors need to include at least one independent dataset to demonstrate that one or multiple parts of their results, e.g., aging-disease associations, GWAS or causal inference results, can be replicated.\nFor the replication data, please consider biobank or atlas data other than UKB, e.g., HuBMAP, All of Us, or MVP. Or alternatively, consider using tissue or disease-specific datasets for replicating specific result of interests, e.g., KPMP for kidney, PsychENCODE for brain, etc.\n(2) Inadequate information disclosed for methods.\nNot enough details were described for key approaches, especially the training of model that imputes age using protein-marker data. Why a local regression between the imputed age and chronological age via lowess function estimation will reveal the true population mean? Details as such, even if referenced from literature, are unclear to most people unless appropriately addressed. Based on the current level of disclosure, it is difficult for readers to understand or trust the methods.\nIn addition, please avoid using computer code-like phrase (i.e., those with an underline, e.g., Line 389-391) in the main text. If the authors feel necessary to include computer codes for the purpose of result reproducibility, please do it within the supplements as a technical protocol.\n(3) Inappropriate arguments.\nThe authors calculated the genetic correlation between the aging statuses of different tissues, and argued that the aging among these tissues was heterogeneous… However, as shown by Figure 3C, heritability of aging in all tissues were low or merely moderate. By my experience, this level of heritability is expected; in this situation, I do not think it is fair to interpret the low genetic correlation here as evidence to the aging heterogeneity, since the heritability itself is low in overall.\nAging is a very complex trait that only a small (or moderate) proportion is genetically determined in each tissue; heterogeneity is true, whereas it is not because of genetic correlation but many contributing environmental factors.\n(4) Exon annotation for prioritizing SNPs from the GWAS result is problematic.\nWhy did the authors use exon SNP annotation as an additional means of prioritizing the significant SNPs in GWAS? To my understanding, the only plausible reason seems to be that exons determine the protein sequence and structure, and thus if any nucleotide substitution in the exons, there will likely to be proteomic changes. However, for purpose of gene prioritization, I do not think the authors should limit the SNPs on exonic regions, although the aging traits were imputed from proteomics. Furthermore, exonic SNPs only constitute 12% of all SNPs in the genome and they are not necessarily non-synonymous.\nTo this end, I would encourage the authors to include intronic SNPs in the vicinity of splice junction sites. These SNPs may alter splicing and thus change the protein structure as well. In addition, please use professional bioinformatics tools or databases to evaluate the regulatory functionality of SNPs rather than using the exon annotation that is too broad and inaccurate.\nHere below I list a few tools or databases for SNP function evaluation\n• HGMD, ClinVar\n• RegulomeDB, HeploReg (TRANSFAC, JASPAR)\n• CADD, DANN, FATHMM-MKL, etc.\n• ExonImpact, regSNPs-splicing, regSNPs-intron\nPlease update this part and also any following parts of the analyses. See if any improvements over the current results and statistical significance.\nIn addition to the major comments above, I also have a number of minor questions as following below. Please also address them.\n(a) Figure 1, e.g., the bottom panel, may not be necessary.\n(b) Figure 2A would be better if changing another way of presenting. Showing only SD increase of delta age is not informative, as people may care more about statistical significance. Consider a volcano plot-like scheme.\n(c) Row labels in Figure 2D have shifted aberrantly.\n(d) Figure 4D is not informative as only showing the location. Try to show direction of effects.\n(e) Had the statistical results of aging status with respect to ABO (Figure 4E-4I) been adjusted for any covariates? If so say so; and if not, please consider all possible covariates and re-perform the analysis.\n(f) Please use uniform scale (e.g., -4 to 6) for the y-axis in Figure 4E-4I.\n(g) Line 280, typo. “Figure 6A” should be “6B”.\n(h) Line 107, inappropriate argument. Rigorously, when a non-significant result, e.g., p-value > 0.05, is derived from some statistics between A and B, we can only say that the current data do not support a significant correlation between them. Whereas this should not be interpreted as A and B have no correlation at all. First, the delta age may not be the best measure of aging; second, the statistical estimation here was based on up to about 50,000 participants, it may be good enough, but we never know before seriously implementing a power analysis. Please keep all discussion scientifically rigorous throughout the paper.\n(i) Line 217-218, inappropriate argument. To my understanding, the genotype-based analysis of ABO is based on the GWAS. Thus, there should not be arguments like “Consistent with the GWAS results…”."}, "author_response": "Response to reviewers’ comments:\nReviewer #1 (Remarks to the Author):\nIn this study, the authors constructed models to predict organ-specific aging using proteomics\ndata, and conducted a thorough investigation searching for genetic determinants of organ-\nspecific aging in large populations. Extensive functional annotations were performed to explore\neffector genes of the identified loci. They also construct polygenic risk scores to assess\nassociations with aging-related diseases for clinical implications. This study is well designed\nand provides very interesting results of genetic architectures for organ aging, which provide\nnew insights into the understanding of the aging process and facilitate the discovery of anti-\naging interventions. The methodology are sound to support the conclusions. The manuscript is\nwell-written. Overall, this is a very nice piece of work, which deserves to be published. I have\nonly a few minor points to potentially improve this work.\n1.Comment: The authors looked up whether the significant loci associated with corresponding\norgan-associated phenotypes in previous GWAS studies. It would be better to check the\nprevious association results for phenotypes related to other organs, since \"no organ system is\nan island\".\nResponse: Thank you for your comment. As you suggested, we checked whether the significant\nloci were associated with phenotypes related to other organs. As shown in current\nSupplementary Table 11, we found that risk SNPs can be annotated not only to the phenotypes\nof the corresponding organs but also to those of other organs, and we provided several examples\nin the Results section (Main text, Line 213-221).\n2.Comment: When comparing with previous GWASs, the authors should present the results of\noverlapping findings.\nResponse: Thank you for your comment. The results of overlapping findings are now shown\nin Supplementary Table 22.\n3.Comment: Methods: The authors should provide the weights of proteins in each model to\nhelp readers understand the importance or contribution of each protein in the model’s\npredictions.\nResponse: Thank you for your comment. We have provided the weights of proteins in\nSupplementary Table 4.\n4.Comment: Specify the method used to assess the relationship between delta ages and\ndiseases, which should be explicitly described in the text.\nResponse: Thank you for your comment. We used logistic regression to assess the relationship\nbetween disease status and delta age. We have provided this information in the Methods section\n(Main Text, Line 468-470).\n5.Comment: The authors should make the GWAS summary data publicly available.\nResponse: Thank you for your comment. We have made the GWAS summary data publicly\navailable in the following link: https://zenodo.org/records/15009392, and this information is\nnow added in the Data availability section (Main text, Line 632-633).\n6.Comment: In line 277, \"supplementary Table 17\" should be corrected to \"supplementary\nTable 16\".\nResponse: Revised as you suggested.\n7.Comment: The X-axis title is missing in Figure 5A.\nResponse: The X-axis title is now added in Figure 5A.\nReviewer #2 (Remarks to the Author):\nZhu and colleagues used Olink plasma proteomics data from the UK Biobank to study the\ngenetic architecture of the multi-organ proteomics age derived from machine learning models.\nProteomics age has been an increasingly used biomarker in aging research. The authors used\nthree different machine-learning models to derive the organ-specific age. Then, they used the\nZ-scores aging biomarkers as phenotypes to run a GWAS and several post-GWAS analyses.\nThe authors made a great effort to analyze the data extensively. The topic is needed in a timely\nmanner and may provide additional insights on previous research endeavors. However, the rigor\nand interpretation of these results need to be further scrutinized. Below are some of my concerns\nand comments.\n1. Comment: Writing: Overall, the writing and flow of the current paper need improvement.\nEspecially for the Abstract, the main messages are obscure, and it seems more like a Method\nsection, not a highlight of the critical findings in this study. For example, the authors state:\n\"Furthermore, we construct polygenic risk scores to assess their association with aging-related\ndiseases. Subsequent Mendelian randomization analysis revealed the potential causal effects of\norgan aging on disease risk status\". I don't think these are the main messages. Instead, they are\nthe specific methods and analyses the authors used in this study. Similarly, for the rest of the\nmanuscript, a clear flow should guide the organization of the writing, including the intuition for\nthe genetic analyses you chose to run.\nResponse: Thank you for the comment! We have carefully revised the manuscript to make sure\nthe writing and flow are clear.\n(1) Specifically, we have highlighted the critical findings of this study in the Abstract as follows\nto be clear:\nHere, using plasma proteomic data from 53,021 participants of European ancestry from the UK\nBiobank, we constructed aging prediction models for 13 organs. GWAS was then performed to\ncharacterize the genetic architecture of the aging for these organs and 125 genomic loci\n(P < 5 × 10–8) were identified. Most of the genetic variants showed organ specificity, while 29\nloci were found to be associated with the aging of multiple organs. Gene prioritization identified\n591 risk genes, which were enriched in biological processes related to organ function. Further\nanalysis showed that the polygenic risk scores of organ aging could be used to predict disease\nstatus. Mendelian randomization analysis revealed potential causality between chronic diseases\n(for example, heart failure and hypertension), modifiable lifestyle factors (smoking and\ndrinking behavior) and organ aging. These findings strengthen our understanding of the genetic\narchitecture of organ aging and illustrate potential strategies for improving human organ health.\n(2) To make sure rest of the manuscript is clear, we have stated the intuition of each analysis in\nthe results section. We have also added a paragraph in line 76-86 to state the organization flow\nas follows to be clear:\nThe overall design of this study is shown in Figure 1. Briefly, using plasma proteomics data\n(2,923 proteins) from 53,021 participants of UKB-PPP, we mapped putative organ-specific\nplasma proteome and used them to train models to predict organ aging for 13 human organs\n(adipose, artery, brain, heart, immune, intestine, kidney, liver, lung, muscle, pancreas, skin)."}, {"round": 2, "reviewer_comments": {"reviewer_2": "The additional comments from the reviewer are attached as a PDF", "reviewer_3": "Overall, I appreciate the authors’ diligent efforts in attempt to resolve the issues in this revision. However, as I evaluated, the issues were partially resolved, not entirely. The remaining issues are still critical. Thus, I do not think this manuscript has been in a good shape for publication unless the remaining issues were seriously taken into consideration.\nHere below I provide my further comments and questions regarding these issues that still remained in the revised manuscript:\n1. “Validation of the aging-disease associations”.\nThe process of model re-training was possibly erroneous. This is because the Elastic net method was involved in the original modeling (using the UKB data). Thus, certain proteins were ought to be selected as features in this step. The correct way of re-training models using other (independent) datasets or cohorts should be fixing on the originally selected features; in specific:\n1) If those selected features also exist in the replication cohort, just use them per se with the exact pre-computed coefficients (in UKB) to impute the aging status (i.e., delta age). In this case, the model does not need and should not be re-trained.\n2) If those selected features largely (but partially) exist in the other cohort, e.g., only a few missing, the prediction model should be re-trained based on the existing features only. This process usually will slightly adjust the coefficient values.\n3) If those features rarely or do not show in the other cohort, you will need to discard this organ (or tissue) aging model as it is technically non-replicable. In this case, you also need to adjust the number of hypotheses accordingly if doing any multiple-hypothesis testing for the modeling results.\nSuch information was not reflected in either the response letter or the currently revised manuscript. To my understanding, the authors re-trained the models de novo in the replication cohorts and thus, this was not a rigorous replication. The results broadly indicated certain organs indeed associate with specific disease (risk), however, nothing to do with the validity of the (original) prediction model established upon the UKB data.\nSuggestion: please follow the standard way of model result replication, as well as provide information about which specific proteins were selected as predictive features with their coefficients.\n2. “Validation of GWAS result”\nIn fact, I did not specifically request that the GWAS signals need to be replicated in my last revision comments. Nonetheless, I do appreciate the authors’ thinking by dividing the large cohort into the northern and southern sub-cohorts, which were with approximately equal number of individuals. This can be a sensible way of finding robust GWAS signals. But since the authors have done this, I have a following question that needs to be clarified:\nI noticed that the author used a lowered threshold of significance (as of p < 0.05) for replication. Why was that?\nTheoretically, since the authors took 10,839 loci to replication, the cutoff of significance level should set to be p < 0.05/10839 (or an FDR corrected by the total number 10,839). Was there any reason for a loosened bar of significance that I might have missed?\n3. “Validation of the causal effects of organ aging on disease status”\nKey information on the associations between SNPs and organ aging (status, delta age) was missing.\nMy understanding is that any causal inference method considers the association between the instrumental variable (IV) and risk factor, together with the association between the IV and disease in certain forms (e.g., summary-statistics). However, how did the authors assess the association of the genetic variants (i.e., IVs) with organ aging (i.e., the risk factor or exposure)?\nIf this was via any method customized or developed by the authors, please provide a description in the manuscript; if via any built-in features in the existing methods, e.g., MR-Egger, please clarify and briefly describe how they did it.\nThe missing of this information raises doubts about the validation (regarding the causal effects). Importantly, again, to qualify a rigorous validation or successful replication, the IVs must be the same.\nDid the authors use exactly the same specific SNPs (of the original MR analysis) in conducting the casual inference analysis in the replication cohorts? If yes, please clarify it in the manuscript; if no, the authors might need to consider re-performing the validation for the causal effects that were derived from the original analysis.\nPlease address all questions above and (further) update the manuscript accordingly. After that, I will re-evaluate whether the answers lead to satisfactory clarification of the issues."}, "author_response": ""}, {"round": 3, "reviewer_comments": {"reviewer_2": "I appreciate the authors’ efforts in revising the manuscript in response to my comments. While they have been generally responsive, I believe some of my concerns remain insufficiently addressed. Below, I outline my perspectives on two specific points that I feel still require further justification:\nComment 2 regarding overfitting: I don't think the author fully understands what an independent test dataset is. Ideally, independent means that the participants were drawn from another study beyond the UK Biobank, but this is not realistic as there is no other biobank that might have this scale. So, what I meant \"indpendent\" test data is that you left out a small proportion of the UK Biobank people at the very beginning (say N=1000, as this sample size is large enough for an objective evaluation of model performance given the large sample size in UK Biobank).\nWhat you have done in this revision is called nested cross-validation, where you have training/validation datasets for model selection, and the test data for within-distribution model evaluation (that is, your cross-vaidated test dataset). However, you might still risk overfitting the cross-validated test data, especially given that you have done so many experiments. The left-held-out \"independent\" test data gives you an unbiased evaluation of your chosen model.\nComment 5 regarding organ-enriched proteins:\nThank you for these comparisons, and I have two follow-up questions:\ni) As you and Oh used the same GTEx data, OLINK proteomics from the UK biobank, and the same criteria, why do you still have some differences for the organ-enriched proteins, say for the brain?\nii) HPA uses various resources, beyond gene expression data only, like GTEx, including proteomics data defined organ-enriched proteins, so in theory, HPA has a better definition of organ-enriched proteins. This is also supported by your comparative results: \"We manually checked all our organ-enriched proteins in HPA, and found that 69.51% (310/446) of them are also supported as enriched in their respective organs in the HPA dataset. 96 of the rest 136 genes were supported as being enhanced (at least four times the mean of all tissue/region/cell types) in their respective organs.\" This suggests that many of the remaining 96 out of 136 genes identified as organ-enriched in your study may not be truly organ-enriched.\nAccording to the Human Protein Atlas (HPA), 'organ-enhanced' proteins are supported by weaker evidence and may not meet the stricter criteria for genuine organ-specific enrichment.\nI don't think GTEx alone should be used to define these proteins as organ-enriched, especially since GTEx data only provides evidence at the gene expression level. More protein-specific resources should be used, like the HPA platform, or more recently developed methodology: https://www.cell.com/cell/fulltext/S0092-8674(25)00286-7", "reviewer_3": "I appreciate the authors’ continuing efforts to address my concerns regarding the scientific quality of results. Based on the latest revised paper, the major issues that I mentioned in the previous revision have been addressed; concerns for scientific clarity, reproducibility, and significance have been alleviated. It is also good to see that authors have provided comprehensive and transparent information with the supplementary data to support the main results.\nHere I have a number of minor suggestions as follows:\n(1) Please think about appropriately highlighting the novelty of this research. Even in this version, this message is insufficiently accredited or discussed. What is the major advantage of this work over others that built or involved proteomics tissue-specific aging models (e.g., Goeminne, L.J.E. et al, 2025, PMID 39488213)? It is the causality study, which presented replicable and significant results. This part is just the limitation of the research of Goeminne, L.J.E. et al. Adding statements of novelty at Introduction and Discussion would be beneficial to publication.\n(2) Although it is a pity that immune system disease were replicated in none of the eight organs (merely nominal p < 0.05), there are five (out of nine) types of diseases for which significant and robust associations were found with respect to seven (out of eight) organs in overall. In my opinion, this result is acceptable because it is understandably the limitation of the follow-up dataset. We could see the method is able to produce robust results, thus lay out a framework for forthcoming studies when the follow-up dataset is enlarged in the future. Nonetheless, one minor issue – regarding these and other statistical results, the authors were saying “multiple testing correction” all the time. Exactly which method? Bonferroni? FDR? This information is in neither Results nor Methods; thus it needs to be clearly mentioned.\n(3) Several supplementary figures are of poor resolution and thus not in a good shape of publication-quality. Especially, Supplementary Figure 8. It is very fuzzy and unable to be viewed clearly. Interestingly, the inset figure of the result replication (i.e., predicted organ-aging with disease) has good visual quality.\n(4) For the MR analysis, it is good that authors added in the missing information that the IVs are the same. However, it still needed to clearly state that the association of the IVs and organ aging was derived from the GWAS results in the paper. I can understand this part because of my study background, but people who are not familiar with MR may not get it at the first sight. Clarifying this point is beneficial to the legibility of the paper.\nI would suggest an acceptance for the present paper, provided that the authors have the manuscript finished in consideration of the issues above."}, "author_response": ""}, {"round": 4, "reviewer_comments": {"reviewer_2": "Thank you to the authors for their response. However, I maintain my position that there is a fundamental misunderstanding regarding nested cross-validation. Based on the explanation provided, it appears that the authors conducted a 5-fold cross-validation that served both for model evaluation and hyperparameter tuning. This setup effectively yields a validation accuracy rather than a proper assessment of generalization performance (cross-validated test and \"independent\" test dataset). Moreover, the authors' referred \"test\" set appears to have been used only once, without being incorporated into a cross-validated framework.\nFor a detailed explanation of nested cross-validation, I recommend referring to this paper by one of the authors of scikit-learn: https://www.sciencedirect.com/science/article/abs/pii/S105381191630595X I believe this point is important, particularly since the authors have gone through multiple revisions, which may have led to some degree of overfitting to the data.\nAlso, I think the authors should make the code publicly available: \"No custom code or mathematical algorithm was used in the methods. All statistical analyses were conducted in Python3 using published libraries and functions.\" This is unacceptable for an open-access journal like Nature Communications for reproducibility.\nThe authors have made a substantial effort in revising the manuscript. I acknowledge that not all published studies adhere to best practices in deriving biological age measures. I will defer to the editor's judgment on whether the manuscript is suitable for publication in its current form."}, "author_response": ""}, {"round": 5, "reviewer_comments": {"reviewer_2": "The authors overall addressed my concern. I checked the code to ensure the nested CV is properly implemented. I found the minor issue:\n- In a strict nested CV, the model chosen in the inner CV (with optimal parameters) is not retrained again on the entire training fold before test evaluation — the inner CV results directly determine the outer test performance.\n- Here, after hyperparameter selection, you refit the model on the whole training fold, which introduces a small optimistic bias.\nI am satisfied with the revision. I do encourage the authors to better understand the best practices in machine learning regarding model generalization and cross-validation for future research."}, "author_response": ""}]

Decision letter