999精品在线视频,手机成人午夜在线视频,久久不卡国产精品无码,中日无码在线观看,成人av手机在线观看,日韩精品亚洲一区中文字幕,亚洲av无码人妻,四虎国产在线观看 ?

一個基于能量的突觸可塑性模型

2022-12-17 03:40:06周一葦陳煥文
生物信息學 2022年4期
關鍵詞:模型

周一葦,陳煥文

(中南大學 自動化學院,長沙 410083)

能量是控制神經元生物物理學機制和神經活動的基礎,神經系統的信息處理離不開大量能量的消耗[1-2]。研究表明靜息狀態下人腦能量消耗占全身的20%,工作狀態下占40%,其中信號傳遞占大腦皮層總能量消耗的80%[3]。這種能量需求大到足以影響人類大腦的設計、功能與進化[4]。在非人類的大腦中能量也能影響大腦功能,例如在小鼠缺氧暴露試驗中就觀察到錐體神經元的突觸阻滯現象與空間記憶缺陷[5]。更進一步地,有研究表明在具有梯度電位和脈沖的神經元中,信噪比、胞體響應速度以及離子通道的數量等方面都與能量有密切關聯[6-7]。提高信噪比需要額外的能量,因為可靠性隨著用于產生信號的隨機事件數而增加。提高帶寬也需要額外的能量,因為增加電導來降低膜時間常數會導致單位時間內需要更強的刺激才能喚起胞體脈沖[6]。改變離子通道數量也需要額外的能量支持,增加鉀離子通道或減少鈉離子通道后,神經元需要消耗更多能量來產生胞體脈沖[7]。由于信噪比、響應速度以及離子通道數量等因素的改變都與電導的變化密切相關,而電導的改變意味著突觸權重發生變化,因此能量也是突觸權重變化的基礎。

能量是影響突觸權重變化的關鍵因素。首先,能量影響突觸前的信息傳遞。能量不僅可以驅動突觸囊泡周期,而且能通過調節突觸囊泡周期的速度來改變囊泡釋放概率,從而影響突觸后權重變化[8-9]。其次,能量約束離子通道動力學以及離子通道數量。由于樹突整合是神經元信息處理過程中主要且昂貴的代謝步驟[10],離子通道的任何變化都可能導致其能量消耗的波動[11],因此能量的上限限制了離子通道的精確類型及其密度[12-13]。綜上所述,能量是突觸權重變化的關鍵因素,給出能量與突觸權重的函數關系是具有價值的。

為系統解釋突觸范圍內能量的變化與突觸可塑性的關系,利用Brian2[14]神經元模擬器構建一個基于能量的突觸可塑性模型。該模型的突觸權重變化方向與幅度取決于正離子進入突觸后膜的瞬時功率。當瞬時功率超過動態功率閾值時突觸權重增長,反之突觸權重下降,增長與下降的幅度由瞬時功率的幅度決定。為檢驗能量模型是否能夠解釋脈沖頻率依賴可塑性以及脈沖時間依賴可塑性這兩種最主要的突觸可塑性,用能量模型復現Sj?str?m等[15]在小鼠視皮層第5層錐體神經元上進行的突觸可塑性實驗。在此基礎上,將能量模型與電壓依賴突觸可塑性模型以及脈沖時刻依賴(Spike-timing dependent plasticity, STDP)突觸可塑性模型進行比較。

1 方 法

1.1 基于能量的突觸可塑性模型建模

使用Python進行網絡仿真,并采用Bono與Clopath基于Brian2[14]開發的詳細生物物理學神經元模型對小鼠視皮層第5層錐體神經元進行模擬[16]。該神經元由1 181個隔室組成,每一個隔室均包含泄露電流、電壓門控型鈉離子通道、延遲整流型鉀離子通道、高閾值鈣離子通道與低閾值鈣離子通道,能夠模擬樹突、軸突和胞體在接受外界刺激時的膜電位變化過程。突觸部分包含α-氨基-3-羥基-5-甲基-4-異惡唑丙酸(α-Amino-3-hydroxy-5-Methyl-4-isoxazole-Propionic Acid, AMPA)受體與N-甲基-D-天冬氨酸(N-methyl-D-aspartic Acid, NMDA)受體,并且具有突觸可塑性,即在外界刺激的情況下,受體的電導可以發生較為持久的變化,也就是突觸權重可以發生持久改變。此處的突觸權重指的是,在突觸受到刺激時突觸后膜受體通道的電導瞬時值與電導最大值的比值。具體公式如下所示:

g(t)=w(t)×gmax

(1)

式(1)中,g(t)是突觸激活時突觸后膜AMPA與NMDA受體的瞬時電導之和,w(t)是突觸權重,gmax是AMPA與NMDA受體的最大電導之和。突觸刺激后AMPA與NMDA受體的電導會依據各自的時間常數指數衰減到零。

在能量模型中突觸權重與能量的關系是,當電流進入突觸后膜的瞬時功率Pm(t)大于可變功率閾值Pth(t)時突觸權重增加,反之突觸權重減少。這是因為在動作電位上升階段,突觸刺激導致AMPA受體、NMDA受體以及其他電壓門控離子通道開放,神經元外部大量的正離子通過這些離子通道順濃度梯度進入突觸后膜,使得突觸后膜的瞬時功率迅速上升,當Pm(t)>Pth(t)時,會導致肌動蛋白細胞骨架和膜動力學改變,使得樹突棘結構和功能的改變,誘發長時程增強(Long Term Potentiation, LTP)[17],突觸權重上升。在動作電位下降階段,由于神經元內的電勢能較高,導致細胞膜內外電勢差較靜息狀態時更小,加之鈉鉀泵消耗能量向細胞膜內外轉運正離子,因此電流涌入突觸后膜的瞬時功率逐漸衰減,當Pm(t)

w(t)=(Pth(t)>θ)×A(t)×Pm(t)

(2)

式中,Pm(t)為突觸后膜的瞬時功率,可變幅值A(t)與功率閾值Pth(t)都是功率Pm(t)的函數(見公式(4)~(8) ),θ為權重變化閾值,取值為0.8 nW,因此突觸權重的變化量w(t)是功率Pm(t)的函數。需要說明的是,突觸權重只有在功率閾值Pth(t)>θ時才發生改變,這是因為只有刺激達到一定水平才能夠激活突觸,并且誘發突觸形態和功能的改變[20-21]。

式中,突觸后膜局部范圍內瞬時功率Pm(t)的計算公式如下:

Pm(t)=Ipost(t)×Vpost(t)×S

(3)

式中,Ipost(t)為突觸后膜單位面積流入的電流,Vpost(t)為突觸后局部膜電位,這兩者皆可由Brian2內置的狀態檢測器直接測出。S為突觸所在隔室的面積,具體參數請查閱參考文獻16。

功率閾值Pth(t)是瞬時功率Pm(t)經過低通濾波后的功率:

(4)

式中,τ為低通濾波的時間常數,其計算公式如下所示:

τ=(36-36×w(t))

(5)

低通濾波的時間常數τ由突觸權重w(t)決定。突觸權重越高,低通濾波時間常數τ越小,功率閾值Pth(t)越接近瞬時功率Pm(t),需要更高頻率或更大幅度的電流刺激才能使瞬時功率Pm(t)超過功率閾值Pth(t)。(5)式中的參數由嘗試法求得。

式(2)中,可變幅值A(t)的公式如下所示:

A(t)=Altp(t)×(Pm(t)-Pth(t))++
Altd(t)×(Pm(t)-Pth(t))-

(6)

式中,Altp(t)與Altd(t)分別為突觸增強與減弱時的可變幅值。公式(Pm(t)-Pth(t))+在Pm(t)>Pth(t)時等于1,反之等于0。(Pm(t)-Pth(t))-則與之相反。突觸權重增強與減弱的幅值不同是符合生物實驗[22]。Altp(t)與Altd(t)都是功率閾值Pth(t)的函數:

Altp(t)=20.7×e-6×Pth(t)+1.4

(7)

Altd(t)=-16×e-6×Pth(t)

(8)

從方程式(6)~(8)中可以看出可變幅值A(t)與功率閾值Pth(t)有關,由于功率閾值Pth(t)是瞬時功率Pm(t)的函數,所以可變幅值A(t)也是瞬時功率Pm(t)的函數。

突觸權重的上升幅度Altp(t)與功率閾值Pth(t)呈指數關系的原因在于,如果突觸后神經元接受單脈沖刺激,由于突觸小泡中釋放的谷氨酸分子數量遠大于突觸后膜受體的數量,因此單脈沖刺激會激活突觸后膜的大部分離子通道[23],導致電流流入神經元的瞬時功率大幅度升高,使得突觸權重的上升幅度Altp(t)隨功率閾值Pth(t)的升高而線性上升。如果突觸后神經元接受的是高頻脈沖刺激,谷氨酸結合位點逐漸接近飽和[24-25],突觸后膜離子通道的電導也接近最大值,此時雖然功率閾值Pth(t)仍在上升,但突觸權重的上升幅度Altp(t)會趨向一個固定值。

突觸權重的下降幅度Altd(t)與功率閾值Pth(t)呈指數關系的原因在于,在突觸刺激剛消失時,突觸后膜電位向靜息電位衰減,由于大部分離子通道仍保持開放狀態,外界的正離子還在通過離子通道進入神經元內部,因此電流流入神經元的瞬時功率Pm(t)僅僅略低于功率閾值Pth(t),且功率閾值Pth(t)保持在較高水平,此時隨著功率閾值Pth(t)的下降,突觸權重的下降幅度Altd(t)仍然接近0。當突觸后膜電位衰減至靜息電位時,突觸后膜的離子通道逐步關閉,此時隨著功率閾值Pth(t)的降低,突觸權重的下降幅度Altd(t)也增大。

1.2 檢驗突觸可塑性的刺激方案

為證明能量是決定突觸可塑性的關鍵因素,構建兩個生物物理學神經元模型(見圖1a),紅點與藍點分別表示基底樹突近端與遠端突觸位點,并仿照Sj?str?m等[15]的生物實驗設計以下兩種刺激方案。

脈沖頻率依賴刺激方案(見圖1b)。為檢驗能量模型能否復現脈沖頻率依賴的突觸可塑性,將兩個神經元的基底樹突近端或遠端通過單個突觸相互連接,分別向兩個神經元的胞體注入3 ms的1 000 pA電流,從而喚起相同頻率的胞體脈沖(0.1、10、20、30、40和50 Hz),每一種頻率下喚起的胞體脈沖總數都為60個。在前—后配對時,突觸前神經元每一個胞體脈沖的起始時間都比突觸后神經元的提前10 ms,后—前配對則相反。突觸權重初始值設為0.5。該刺激方案分別在樹突近端與遠端重復10次,每一次均重新選取突觸位置。利用Brian2內含的狀態檢測器來檢測突觸權重與突觸后膜瞬時功率的變化,得到圖2a~2f。圖2a與2d的陰影部分表示一個標準差的范圍。

脈沖時間依賴刺激方案(見圖1c)。為檢驗能量模型能否復現脈沖時間依賴的突觸可塑性,將兩個神經元的基底樹突近端或遠端通過單個突觸相互連接,通過向胞體注入3 ms的1 000 pA電流的方式,在兩個神經元上都喚起頻率為20 Hz的胞體脈沖,2個神經元胞體脈沖起始時間之間具有時間差(1、2.5、5、7.5、10、20和30 ms),并統計60次脈沖后突觸權重的變化。前—后配對時突觸前神經元先釋放一個胞體脈沖,在經過間隔時間(1、2.5、5、7.5、10、20和30 ms)后突觸后神經元釋放一個胞體脈沖。后—前配對則相反。突觸初始權重為0.5。該刺激方案分別在樹突近端與遠端重復10次,每次都重新選取突觸位置。用狀態監測器記錄突觸權重變化與突觸后膜瞬時功率的變化,得到圖2g~2k。圖2g的陰影部分表示一個標準差的范圍。

圖1 神經元連接及刺激方案示意圖Fig.1 Schematic diagram of neuron connection and stimulation scheme

2 結 果

2.1 突觸可塑性檢驗

為研究能量模型是否能復現脈沖頻率依賴可塑性,將兩個生物物理學神經元相連接,采用Sj?str?m等[15]設計的脈沖頻率依賴刺激方案。圖2a展示了在前后-配對的脈沖頻率刺激方案下,突觸權重隨配對頻率的變化。紅線與藍線分別是樹突近端與遠端的突觸響應平均值,陰影部分代表一個標準差的范圍。誤差棒表示Sj?str?m等[15]通過生物實驗得到的誤差范圍,黑色實線是所有突觸的響應平均值,虛線是初始權重0.5??梢钥闯觯谇啊笈鋵r,無論突觸位于基底樹突近端還是遠端,刺激頻率越大,突觸前神經元連接至突觸后神經元的突觸權重的上升幅度也越大,并且刺激后突觸權重的平均值都落在生物實驗的誤差范圍內。這是因為在突觸刺激與突觸后神經元的胞體脈沖共同作用下,突觸后神經元的離子通道打開,大量正離子順離子濃度梯度進入突觸后膜,導致注入突觸后膜的電流增大,由于突觸后膜的瞬時功率等于進入突觸后膜的瞬時電流乘以細胞膜的電位差,因此突觸后膜的瞬時功率上升。當瞬時功率Pm(t)高于功率閾值Pth(t)時,突觸上的肌動蛋白細胞骨架和膜動力學發生改變,突觸權重上升[17]。

圖2b與2c進一步解釋了在脈沖頻率依賴刺激方案中以30 Hz進行前—后配對時,突觸后膜的瞬時功率Pm(t)與突觸權重w(t)隨時間變化的關系。圖2b的綠線是突觸后膜的瞬時功率Pm(t),黑線是功率閾值Pth(t)。由圖2b與2c可知,在0至60 ms時,由于沒有電流刺激,功率閾值Pth(t)低于權重變化閾值θ,突觸權重不發生改變。在60 ms至70 ms階段,在突觸刺激與突觸后神經元胞體脈沖共同作用下,突觸后膜的瞬時功率Pm(t)高于功率閾值Pth(t),導致突觸權重上升。在70 ms~110 ms階段,由于離子通道逐漸關閉,Pm(t)

脈沖頻率依賴刺激方案的后—前配對結果如圖2d所示。由圖可見如果兩個神經元的胞體脈沖頻率相同且突觸前神經元的胞體脈沖比突觸后神經元落后10 ms,無論突觸位于基底樹突的近端或遠端,突觸前神經元連接至突觸后神經元的突觸權重在刺激頻率小于35 Hz時都會下降,只有在刺激頻率大于35 Hz時才會逐步上升。這是因為突觸后神經元在胞體脈沖后進入去火期,在此期間進行突觸刺激只能打開少許突觸后膜的離子通道,引起少量正離子流入,即突觸后膜的瞬時功率Pm(t)低于功率閾值Pth(t),所以無法引起突觸權重的改變甚至會導致誘導LTD。

圖2e與2f進一步展示了在脈沖頻率依賴刺激方案中以30 Hz進行后—前配對時,突觸后膜的瞬時功率Pm(t)與突觸權重w(t)變化的關系。在60 ms之前,由于沒有突觸刺激與胞體脈沖,突觸后膜的瞬時功率Pm(t)為零,突觸權重無變化。在60 ms時,突觸后神經元產生胞體脈沖,但由于樹突棘頸部的電阻很大,削弱了胞體脈沖的影響,因此突觸后膜只有極少量離子通道打開,瞬時功率Pm(t)與功率閾值Pth(t)幾乎沒變化,突觸權重也不改變。在70 ms時,突觸刺激引起突觸后膜大量離子通道打開,瞬時功率Pm(t)高于功率閾值Pth(t),突觸權重上升,但由于沒有誘發突觸后神經元的胞體脈沖,在突觸刺激撤去后突觸權重持續下降,在重復5次配對刺激后,突觸權重從初始值0.5下降到0.495。

其次,為研究能量模型是否能復現脈沖時間依賴可塑性,對能量模型采用Sj?str?m等[15]設計的脈沖時間依賴刺激方案。由圖2g可見,只有在采取前—后配對且脈沖時間間隔在10 ms左右時突觸權重增加,在后—前配對時突觸權重的變化量為負值。隨著間隔時間的增大,無論是前—后配對還是后—前配對,突觸權重的變化量都趨于0。這是因為在后—前配對時,突觸后神經元的胞體脈沖始終先于突觸刺激,由于胞體脈沖持續時間較短(1~2 ms)[26],因此兩者無法相互疊加,導致突觸后膜的離子通道只有少量開放,瞬時功率Pm(t)低于功率閾值Pth(t),導致LTD。在前—后配對時,由于突觸刺激持續時間較長,突觸后神經元的胞體脈沖與突觸刺激能夠同時作用在離子通道上,使得大量正離子進入突觸后膜,瞬時功率Pm(t)高于功率閾值Pth(t),導致LTP。當間隔時間超過20 ms后,無論是前—后配對還是后—前配對,突觸刺激和胞體脈沖對離子通道的作用都無法疊加,瞬時功率Pm(t)始終低于功率閾值Pth(t)且接近于0,因此突觸權重的改變量為負值并且趨向0。

圖2h~2k分別展示了在脈沖時間依賴刺激方案中胞體脈沖間隔時間為+10 ms與-10 ms時,突觸后膜的瞬時功率Pm(t)與突觸權重w(t)隨時間的變化關系。由圖2h與2i可知,當脈沖間隔時間為+10 ms時,在60 ms~70 ms,突觸后神經元產生胞體脈沖與突觸刺激共同誘導大量離子通道開放,瞬時功率Pm(t)高于功率閾值Pth(t),重復5次配對后突觸權重上升到0.509。圖2j與2k則顯示當脈沖間隔時間為-10 ms時,在60 ms突觸后神經元產生胞體脈沖,但由于沒有與突觸刺激相配對,因此瞬時功率Pm(t)很小,權重無變化。在70 ms時突觸刺激使瞬時功率Pm(t)上升,由于沒喚起突觸后神經元的動作電位,因此突觸刺激消失后Pm(t)長時間低于功率閾值Pth(t),重復刺激5次后突觸權重最終下降到0.488。

圖2 脈沖頻率以及脈沖時間依賴可塑性檢驗Fig.2 Pulse frequency and pulse time dependent plasticity test

2.2 模型比較

為了證明基于能量的突觸可塑性模型與其他模型相比具有優越性,比較能量模型、電壓依賴可塑性模型[16]以及STDP模型[27]在脈沖頻率依賴可塑性實驗和脈沖時間依賴可塑性實驗中的響應,以Sj?str?m等生物實驗的誤差范圍[15]作為評判模型優越性的標準。

圖3a與3b顯示了三種突觸可塑性模型在脈沖頻率依賴可塑性實驗下的響應,其中紅線、黃線與藍線分別表示能量模型、電壓模型與STDP模型在基底樹突近端與遠端的突觸響應平均值,誤差棒表示Sj?str?m等[15]通過生物實驗得到的誤差范圍,黑色虛線是初始權重0.5。由圖3a可知,在前—后配對情況下,能量模型與電壓模型的響應僅在刺激頻率接近0 Hz時略高于生物實驗的平均值,在其余刺激頻率下都與生物實驗的平均值相似。相比之下,STDP模型的響應只在刺激頻率大于20 Hz時落在誤差范圍內,在低頻刺激下都落在誤差范圍外。再由圖3b可知,在后—前配對的情況下,電壓模型與STDP模型的響應在低頻段(≤20Hz)都貼合實驗值,在高頻段(>20Hz)兩個模型的響應遠離實驗數據。與之相比,能量模型的響應在刺激頻率為0 Hz時離誤差值較遠,但在高頻段更接近生物實驗的誤差范圍。圖3a與3b共同證明了在脈沖頻率依賴可塑性的實驗條件下,能量模型相較公認的電壓依賴可塑性模型以及傳統的STDP模型更具有優越性。

其次,對比脈沖時間依賴可塑性實驗的結果(見圖3c),紅線、黃線、藍線以及誤差棒的意義與圖3a與2b相同。由圖可見,能量模型的響應僅在間隔時間接近0 ms時低于實驗值,在間隔時間大于或小于0 ms時都處在誤差范圍內。而電壓模型和STDP模型的響應只在間隔時間等于-10 ms或10 ms時貼合實驗值,其余部分均不在誤差范圍內。因此證明在脈沖時間依賴可塑性的實驗條件下,能量模型同樣具有優越性。

圖3 模型對比分析Fig.3 Model comparison analysis

3 結 論

與其他的突觸可塑性模型相比,本文提出的基于能量的突觸可塑性模型有兩個主要優勢。第一點是明確了能量是影響突觸可塑性的關鍵因素。影響功能連接模式的因素眾多,為解釋生物實驗中觀察到的突觸可塑性,部分模型將多種因素納入考慮[16,27-29],例如STDP模型[27]就需要考慮突觸前和突觸后兩個神經元的脈沖間隔。相比之下,我們構建的突觸可塑性模型僅通過對突觸后范圍內能量變化的分析,就能夠再現脈沖頻率依賴與脈沖時間依賴的突觸可塑性,并且比電壓依賴的突觸可塑性模型[16]和經典的STDP模型[27]更貼合真實生物實驗的結果。第二點是能量模型能夠從生物物理學角度定量解釋突觸權重變化的過程。目前較為流行的STDP模型[27]、Hebbian模型[28]以及BCM模型[29]都只是描述神經元脈沖響應這一現象與突觸權重之間的關聯,并不涉及神經元內在的變化。而能量模型則是根據突觸范圍內的能量變化解釋突觸權重變化的過程,具有生物物理意義。

能量模型證明能量是解釋突觸可塑性的關鍵所在,這加深了我們對突觸可塑性機制的理解,為學習和記憶的進一步研究提供了一個新的思路。接下來還可以從兩方面對能量模型進行改進。首先,該模型雖然通過采用可變功率閾值Pth(t)的方式來避免硬邊界設計,但在持續接受不符合生物實際的高頻刺激時,突觸權重有可能超出正常范圍。第二,該模型的突觸權重在瞬時功率與功率閾值都回歸靜息狀態后就維持在一個定值。而在真實的神經元中,突觸長時間未接受刺激時,其突觸權重應該隨時間逐漸衰減。可以考慮添加其他機制來模擬這個過程。

猜你喜歡
模型
一半模型
一種去中心化的域名服務本地化模型
適用于BDS-3 PPP的隨機模型
提煉模型 突破難點
函數模型及應用
p150Glued在帕金森病模型中的表達及分布
函數模型及應用
重要模型『一線三等角』
重尾非線性自回歸模型自加權M-估計的漸近分布
3D打印中的模型分割與打包
主站蜘蛛池模板: 午夜国产大片免费观看| 国产精品入口麻豆| 97视频在线观看免费视频| 992tv国产人成在线观看| 亚洲国产高清精品线久久| 国产丰满成熟女性性满足视频 | 18黑白丝水手服自慰喷水网站| 美女啪啪无遮挡| 中文国产成人精品久久一| 一区二区在线视频免费观看| 8090午夜无码专区| 日韩精品一区二区深田咏美| 99久久99这里只有免费的精品| 成人午夜福利视频| 狠狠色综合久久狠狠色综合| 亚洲一区二区成人| 国产亚洲精品自在线| 丁香婷婷激情网| 国产精品福利在线观看无码卡| 91在线精品免费免费播放| 国产农村妇女精品一二区| 国产区在线看| 国产成本人片免费a∨短片| 伊大人香蕉久久网欧美| 久久国产精品夜色| 国产微拍精品| 国产一区二区丝袜高跟鞋| а∨天堂一区中文字幕| 国产91无码福利在线| 亚洲国产日韩在线观看| 亚洲有无码中文网| 综合网天天| 欧美成人综合视频| 欧美一级夜夜爽| 国产一级裸网站| 在线网站18禁| 思思99思思久久最新精品| 一区二区欧美日韩高清免费| 国产导航在线| 国模极品一区二区三区| 中字无码av在线电影| 青青草原国产| 2021国产精品自拍| 97久久精品人人| 国产美女一级毛片| 日韩人妻无码制服丝袜视频| 国产流白浆视频| 国产一区二区三区免费观看| 国产精品护士| 国产精品流白浆在线观看| 国产精品手机在线播放| 亚洲精品无码专区在线观看| 片在线无码观看| 性视频一区| 精品久久蜜桃| 亚洲一区免费看| 久久这里只有精品23| 重口调教一区二区视频| 国产高清毛片| 国产欧美日韩在线在线不卡视频| 无码专区国产精品一区| 国产真实自在自线免费精品| 国产成人免费视频精品一区二区| 91原创视频在线| 成人综合网址| 久久99国产精品成人欧美| 熟女成人国产精品视频| 伊人查蕉在线观看国产精品| 国产第三区| 三级毛片在线播放| 性色一区| 国产专区综合另类日韩一区| 国产乱码精品一区二区三区中文| 色综合国产| 沈阳少妇高潮在线| h视频在线播放| 午夜精品一区二区蜜桃| 91人人妻人人做人人爽男同| 久久特级毛片| 狠狠色噜噜狠狠狠狠奇米777| 四虎国产在线观看| 91九色国产porny|