ループとROC曲線作成の練習②後編（ラベル・値ラベルの抽出、結果のエクセルへの貼付け）

このブログでは、統計解析ソフトStataのプログラミングのTipsや便利コマンドを紹介しています．

Facebook groupでは、ちょっとした疑問や気づいたことなどを共有して貰うフォーラムになっています．ブログと合わせて個人の学習に役立てて貰えれば幸いです．

さて、今回は、前の記事の続きです．今回はサブグループの結果を取り出す作業を一緒にやっていってみましょう．

3．サブグループの結果を取り出してみる
4．エクセルへの貼り付け→表の完成
まとめ

3．サブグループの結果を取り出してみる

さて、同じサンプルデータを使って、年齢と性別のサブグループ別のmarker_aとmarker_bのカットオフ、感度、特異度を取り出してみましょう．

前回とは違って、年齢のカテゴリーの値ラベルを作ってみます．このとき、文字列にスペースを入れないようにするのですが、これは後から意味が分かります．

import delimited "C:\...\sampledata_cutpt.csv", clear  /* ...のところは作業directoryを指定してください */
recode age (min/55=1) (55/65=2) (65/max=3), gen(agecat)
label define agecat 1 "Age<55" 2 "Age55-65" 3 "Age>=65"
label values agecat agecat
encode sex, gen(sex_c)

さて、ここから核になるのは、以下のsyntaxです．outcomeに対するmarker_aの分類性能をみたいので、Youden indexで求めたカットオフ値における感度と特異度を取り出して表示します．

cutptというコマンドによって行った演算の結果を取り出すには以下のところをみて、back tickとアポストロフィで括ってあげればOKです．

cutpt outcome marker_a, noadjust youden
local sens = `e(sens)'*100
local spec = `e(spec)'*100
disp "marker_a", %3.1f `e(cutpoint)', %3.0f `sens`v'`m'', %3.0f `spec`v'`m''

Youden indexにより求めた最適カットオフ値は16であり、感度は61％、特異度は66％とわかる

これを性別（sex_c: 1/2）と年齢カテゴリ（agecat: 1/2/3）別に表示させて最後は表にまとめたいですよね．そこで私はいつも、ここでスペース区切りで結果を表示させて最後にエクセルにコピーして、「区切りの位置」で整える、という技を使って結果を楽に取り出しています．

forvalues m = 1/2 {
	cutpt outcome marker_a if sex_c==`m', noadjust youden
	local sens`m' = `e(sens)'*100
	local spec`m' = `e(spec)'*100
	disp "marker_a",`m', "total", %3.1f `e(cutpoint)', %3.0f `sens`m'', %3.0f `spec`m'' 
   }

このようにすると、すべての年齢におけるmarker_aの男女別の効果を現すことができます．これにさらに年齢層別の結果を後ろにくっつけたい場合、

forvalues m = 1/2 {
	cutpt outcome marker_a if sex_c==`m', noadjust youden
	local sens`m' = `e(sens)'*100
	local spec`m' = `e(spec)'*100
	disp "marker_a",`m', "total", %3.1f `e(cutpoint)', %3.0f `sens`m'', %3.0f `spec`m'' 
  forvalues n = 1/3 {
	cutpt outcome marker_a if sex_c==`m' & agecat==`n', noadjust youden
	local sens`m'`n' = `e(sens)'*100
	local spec`m'`n' = `e(spec)'*100
	disp "",`m', `n' , %3.1f `e(cutpoint)', %3.0f `sens`m'`n'', %3.0f `spec`m'`n''
     }
  }

これだと、途中のコマンドが入ってきて綺麗にコピーできないので、他をすべて表示させないように”quietly”で全体を覆って、結果だけnoisに表示させましょう．

quietly {
forvalues m = 1/2 {
	cutpt outcome marker_a if sex_c==`m', noadjust youden
	local sens`m' = `e(sens)'*100
	local spec`m' = `e(spec)'*100
	nois disp "marker_a",`m', "total", %3.1f `e(cutpoint)', %3.0f `sens`m'', %3.0f `spec`m'' 
  forvalues n = 1/3 {
	cutpt outcome marker_a if sex_c==`m' & agecat==`n', noadjust youden
	local sens`m'`n' = `e(sens)'*100
	local spec`m'`n' = `e(spec)'*100
	nois disp "",`m', `n' , %3.1f `e(cutpoint)', %3.0f `sens`m'`n'', %3.0f `spec`m'`n''
     }
  }   
}

これをmarker_a, marker_bともに行いたいので、

quietly {
  foreach v of varlist marker_a marker_b {	
	forvalues m = 1/2 {
    	  cutpt outcome `v' if sex_c==`m', noadjust youden
	  local sens`v'`m' = `e(sens)'*100
	  local spec`v'`m' = `e(spec)'*100
	  noi disp "`v'",`m', "total", %3.1f `e(cutpoint)', %3.0f `sens`v'`m'', %3.0f `spec`v'`m'' 
		forvalues n = 1/3 {
	 	  cutpt outcome `v' if sex_c==`m' & agecat==`n', noadjust youden
		  local sens`v'`m'`n' = `e(sens)'*100
		  local spec`v'`m'`n' = `e(spec)'*100
		  noi disp "`v'",`m', `n' , %3.1f `e(cutpoint)', %3.0f `sens`v'`m'`n'', %3.0f `spec`v'`m'`n''
		}
	}
  }
}

さらに欲張って、値ラベルを表示させてみます．今のままだと、agecatが1,2,3としか表示されませんので、値ラベルで表示させるために、displayさせるときに、

“`: label (varname) #’ “

（）の中には変数名、そのうちの何番目の変数かを＃が表します．なのでたいていは＃のところやvarnameのところにマクロが入ってきます．（その場合には`# ‘で表現します．ラベルの何番目、と数字で指定する場合には`’入れる必要はないですのでご注意ください．）

qui{
	foreach v of varlist marker_a marker_b {	
		forvalues m = 1/2 {
			cutpt outcome `v' if sex_c==`m', noadjust youden
			local sens`v'`m' = `e(sens)'*100
			local spec`v'`m' = `e(spec)'*100
			noi disp "`v'","`: label (sex_c) `m''", "total", %3.1f `e(cutpoint)', %3.0f `sens`v'`m'', %3.0f `spec`v'`m'' 
			forvalues n = 1/3 {
				cutpt outcome `v' if sex_c==`m' & agecat==`n', noadjust youden
				local sens`v'`m'`n' = `e(sens)'*100
				local spec`v'`m'`n' = `e(spec)'*100
				noi disp "`v'","`: label (sex_c) `m''", "`: label (agecat) `n''" , %3.1f `e(cutpoint)', %3.0f `sens`v'`m'`n'', %3.0f `spec`v'`m'`n''
			}
		}
	}
}

ここまでのプロセスが普通にできるようになるとかなりStataのプログラミングに自信が持てると思います．

もう一度復習ですが、値ラベルを取り出したいときはこれ！

“`: label (varname) # ‘ “

ちなみに変数ラベルを取り出したいときは、

“`:variable label varname ‘ “