12

プロポーションのオブジェクトから信頼区間を作成する既存の関数はありますか（私の場合、パッケージsvyby内のバイナリアイテムのクロス集計）。私はグループ間の比率を比較することがよくありますが、信頼区間を抽出できる関数（ではなくsurvey調査関数を使用）があると非常に便利です。以下の例は、私が達成したいことを示しています。svycipropconfint

データを読み込む

library(survey)
library(weights)
data(api)
apiclus1$both<-dummify(apiclus1$both)[,1]#Create dummy variable
dclus1<-svydesign(id=~dnum, weights=~pw, data=apiclus1, fpc=~fpc)

stype全体で変数「両方」の比率を比較するsvybyオブジェクトを作成します

b<-svyby(~both, ~stype, dclus1, svymean)
confint(b)#This works, but svyciprop is best in  other cases, especially when proportion is close to 0 or 1
svyciprop(b)#This requires that you specify each level and a design object

関数を作成することは可能でしょうか（たとえば、を使用byCI(b,method="likelihood")するのと同じことを実現しますか？基本的にオブジェクトの各レベルを通過して信頼区間を作成する必要があります。これまでの試みは失敗しました。confint(b)svycipropsvyby

svyby()これを回避する別の方法があるかもしれませんが、それは迅速で直感的であるため、私は使用するのが好きです。

4

2 に答える 2

12

svyby()にはvartype=、サンプリングの不確実性をどのように指定するかを指定する引数があります。vartype="ci"信頼区間を取得するために使用します。

svyby(~I(ell>0),~stype,design=dclus1, svyciprop,vartype="ci",method="beta")

これが各レベルを手動で行うのと同じになることを確認するのは簡単です。

confint(svyciprop(~I(ell>0), design=subset(dclus1,stype=="E"),method="beta"))
于 2013-01-02T21:19:14.807 に答える
2

興味深い..これら2つのコマンドは同じ結果を与えるべきではありません..最初のコマンドはおそらくエラーまたは警告をスローするはずです:

svyby( ~both , ~stype , dclus1 , svyciprop , method = 'likelihood' )
svyby( ~both , ~stype , dclus1 , svymean )

この問題についてラムリー博士に警告することをお勧めします - の 80 行目付近のコードは、おそらく内部で動作surveyby.Rするようにわずかに変更される可能性があります..これについて彼に連絡する前に、すべてを注意深く読んでくださいsvycipropsvyby

とにかく、これはあなたの問題を解決するかもしれない一時的な解決策です

# create a svyby-like function specific for svyciprop
svyciby <- 
    function( formula , by , design , method = 'likelihood' , df = degf( design ) ){

        # steal a bunch of code from the survey package's source
        # stored in surveyby.R..
        byfactors <- model.frame( by , model.frame( design ) , na.action = na.pass )
        byfactor <- do.call( "interaction" , byfactors )
        uniquelevels <- sort( unique( byfactor ) )
        uniques <- match( uniquelevels , byfactor )
        # note: this may not work for all types..
        # i only tested it out on your example.

        # run the svyciprop() function on every unique combo
        all.cis <-
            lapply( 
                uniques , 
                function( i ){

                    svyciprop( 
                        formula , 
                        design[ byfactor %in% byfactor[i] ] ,
                        method = method ,
                        df = df
                    )
                }
            )

        # transpose the svyciprop confidence intervals
        t.cis <- t( sapply( all.cis , attr , "ci" ) )

        # tack on the names
        dimnames( t.cis )[[1]] <- as.character( sort( unique( byfactor ) ) )

        # return the results
        t.cis
    }

# test out the results
svyciby( ~both , ~stype , dclus1 , method = 'likelihood' )
# pretty close to your b, but not exact (as expected)
confint(b)
# and this one does match (as it should)
svyciby( ~both , ~stype , dclus1 , method = 'mean' , df = Inf )
于 2013-01-01T18:51:35.340 に答える