Load libraries.

library(dagitty)
library(ggdag)

Attaching package: ‘ggdag’

The following object is masked from ‘package:stats’:

    filter
library(dosearch)

Create a DAG of the experiments. The frame is ignored for simplicity. The DAG applies to each frame set. The following notations are used:

```r
g <- dagify(Y ~ D + R,
            R ~ C,
            C ~ D,
            exposure = \D\,
            outcome = \Y\,
            coords = list(x = c(D = 1, C = 2, R = 3, Y = 4),
                y = c(D = 1, C = 2, R = 1.5, Y = 1))
            )

ggdag(g)

<!-- rnb-source-end -->

<!-- rnb-chunk-end -->


<!-- rnb-text-begin -->


Next, let's see how to do the front-door adjustment for the indirect path. 


<!-- rnb-text-end -->


<!-- rnb-chunk-begin -->


<!-- rnb-source-begin eyJkYXRhIjoiYGBgclxuZGF0YTEgPC0gXCJQKEMsIFIsIFkpXCJcblxucXVlcnkxID0gXCJQKFkgfCBkbyhDKSlcIlxuXG5ncmFwaDEgPSBcIlxuQyAtPiBSXG5SIC0+IFlcblUgLT4gWVxuVSAtPiBDXG5cIlxuXG4jIGNvbXB1dGVcbmZyb250ZG9vciA9IGRvc2VhcmNoKGRhdGExLCBxdWVyeTEsIGdyYXBoMSlcblxuIyBjb252ZXJ0IHRvIFJtYXJrZG93biBlcXVhdGlvblxuY2F0KHBhc3RlKFwiJCRcIiwgZnJvbnRkb29yJGZvcm11bGEsIFwiJCRcIikpXG5gYGAifQ== -->

```r
data1 <- "P(C, R, Y)"

query1 = "P(Y | do(C))"

graph1 = "
C -> R
R -> Y
U -> Y
U -> C
"

# compute
frontdoor = dosearch(data1, query1, graph1)

# convert to Rmarkdown equation
cat(paste("$$", frontdoor$formula, "$$"))
$$ \sum_{R}\left(p(R|C)\sum_{C}\left(p(C)p(Y|R,C)\right)\right) $$

The adjustment equation is: \[ \sum_{R}\left(p(R|C)\sum_{C^{\prime}}\left(p(C^{\prime})p(Y|R,C^{\prime})\right)\right) \]

We apply this equation to both the treatment and control group. In the computation, we can focus on four groups: 2 (treatment vs. control) by 2 (2 frame sets). I will use the treatment group as an example. The control group follows the same logic.

\[\begin{eqnarray} p(F=f|D=1)p(C=c|F=f,D=1)\left[\sum_{R}\left(p(R|C)\sum_{C^{\prime}}\left(p(C^{\prime})p(Y|R,C^{\prime})\right)\right)\right]_{D=1\&F=f} \end{eqnarray} \]

The subscript \(D=1\) means we evaluate the term for the treatment group. For the treatment group, we can calculate the inner term for each frame set and then combine their effects. The term in the bracket needs special attention as \(\sum_{C^{\prime}}\left(p(C^{\prime})p(Y|R,C^{\prime})\right)\) is evaluated for both choices when at different \(R\).

The procedure for evaluation is like this. For each frame set group, we the following:

  1. Calculate this term \(\sum_{C^{\prime}}\left(p(C^{\prime})p(Y|R,C^{\prime})\right)\), at \(R=1\) and \(R=0\). This gives use \(P(Y|do(R)=1)\) and \(P(Y|do(R)=0)\).
  2. With the values, we calculate the value of \(P(Y|do(C)=1)\) and \(P(Y|do(C)=2)\) as: \[\begin{eqnarray*} P(Y|do(C)=1) = \sum_{R}p(R|C=1)P(Y|do(R^{\prime})=R)\\ P(Y|do(C)=2) = \sum_{R}p(R|C=2)P(Y|do(R^{\prime})=R) \end{eqnarray*} \]

With values of \(P(Y|do(C))\), we combine them to get the total path effect at \(F=f\) (I ignore the \(do(F)=f\) as \(F\) is randomized). \[P(Y|do(C),F=f)= \sum_{c}p(C=c|F=f)P(Y|do(C)=c)\]

Given the per-frame set group values, we combine to get the total indirect effect for the treatment group: \[P(Y|do(D)=1,do(C))=\sum_{f}P(F=f|D=1)P(Y|do(C),F=f)\]

LS0tDQp0aXRsZTogIkRhdGEgQW5hbHlzaXMgZm9yIEFSIEV4cGVyaW1lbnQiDQpvdXRwdXQ6IGh0bWxfbm90ZWJvb2sNCi0tLQ0KDQpMb2FkIGxpYnJhcmllcy4gDQpgYGB7ciwgd2FybmluZz1GQUxTRX0NCmxpYnJhcnkoZGFnaXR0eSkNCmxpYnJhcnkoZ2dkYWcpDQpsaWJyYXJ5KGRvc2VhcmNoKQ0KYGBgDQoNCkNyZWF0ZSBhIERBRyBvZiB0aGUgZXhwZXJpbWVudHMuIFRoZSBmcmFtZSBpcyBpZ25vcmVkIGZvciBzaW1wbGljaXR5LiBUaGUgREFHIGFwcGxpZXMgdG8gZWFjaCBmcmFtZSBzZXQuIFRoZSBmb2xsb3dpbmcgbm90YXRpb25zIGFyZSB1c2VkOiANCg0KKiBEIC0gQVIgdnMuIGltYWdlcyAodHJlYXRtZW50KS4NCiogQyAtIENob2ljZSBvZiBnbGFzc2VzLiANCiogUiAtIFJpZ2h0IGZyYW1lIG9yIHdyb25nIG9uZSwgZ2l2ZW4gdGhlIGNob2ljZS4gDQoqIFkgLSBQcm9kdWN0IHJldHVybnMuIA0KYGBge3J9DQpnIDwtIGRhZ2lmeShZIH4gRCArIFIsDQogICAgICAgICAgICBSIH4gQywNCiAgICAgICAgICAgIEMgfiBELA0KICAgICAgICAgICAgZXhwb3N1cmUgPSAiRCIsDQogICAgICAgICAgICBvdXRjb21lID0gIlkiLA0KICAgICAgICAgICAgY29vcmRzID0gbGlzdCh4ID0gYyhEID0gMSwgQyA9IDIsIFIgPSAzLCBZID0gNCksDQogICAgICAgICAgICAgICAgeSA9IGMoRCA9IDEsIEMgPSAyLCBSID0gMS41LCBZID0gMSkpDQogICAgICAgICAgICApDQoNCmdnZGFnKGcpDQoNCmBgYA0KDQpOZXh0LCBsZXQncyBzZWUgaG93IHRvIGRvIHRoZSBmcm9udC1kb29yIGFkanVzdG1lbnQgZm9yIHRoZSBpbmRpcmVjdCBwYXRoLiANCg0KYGBge3J9DQpkYXRhMSA8LSAiUChDLCBSLCBZKSINCg0KcXVlcnkxID0gIlAoWSB8IGRvKEMpKSINCg0KZ3JhcGgxID0gIg0KQyAtPiBSDQpSIC0+IFkNClUgLT4gWQ0KVSAtPiBDDQoiDQoNCiMgY29tcHV0ZQ0KZnJvbnRkb29yID0gZG9zZWFyY2goZGF0YTEsIHF1ZXJ5MSwgZ3JhcGgxKQ0KDQojIGNvbnZlcnQgdG8gUm1hcmtkb3duIGVxdWF0aW9uDQpjYXQocGFzdGUoIiQkIiwgZnJvbnRkb29yJGZvcm11bGEsICIkJCIpKQ0KYGBgDQpUaGUgYWRqdXN0bWVudCBlcXVhdGlvbiBpczogDQokJCBcc3VtX3tSfVxsZWZ0KHAoUnxDKVxzdW1fe0Nee1xwcmltZX19XGxlZnQocChDXntccHJpbWV9KXAoWXxSLENee1xwcmltZX0pXHJpZ2h0KVxyaWdodCkgJCQNCg0KV2UgYXBwbHkgdGhpcyBlcXVhdGlvbiB0byBib3RoIHRoZSB0cmVhdG1lbnQgYW5kIGNvbnRyb2wgZ3JvdXAuIEluIHRoZSBjb21wdXRhdGlvbiwgd2UgY2FuIGZvY3VzIG9uIGZvdXIgZ3JvdXBzOiAyICh0cmVhdG1lbnQgdnMuIGNvbnRyb2wpIGJ5IDIgKDIgZnJhbWUgc2V0cykuIEkgd2lsbCB1c2UgdGhlIHRyZWF0bWVudCBncm91cCBhcyBhbiBleGFtcGxlLiBUaGUgY29udHJvbCBncm91cCBmb2xsb3dzIHRoZSBzYW1lIGxvZ2ljLiANCg0KJCRcYmVnaW57ZXFuYXJyYXl9DQpwKEY9ZnxEPTEpcChDPWN8Rj1mLEQ9MSlcbGVmdFtcc3VtX3tSfVxsZWZ0KHAoUnxDKVxzdW1fe0Nee1xwcmltZX19XGxlZnQocChDXntccHJpbWV9KXAoWXxSLENee1xwcmltZX0pXHJpZ2h0KVxyaWdodClccmlnaHRdX3tEPTFcJkY9Zn0NClxlbmR7ZXFuYXJyYXl9DQokJA0KDQpUaGUgc3Vic2NyaXB0ICREPTEkIG1lYW5zIHdlIGV2YWx1YXRlIHRoZSB0ZXJtIGZvciB0aGUgdHJlYXRtZW50IGdyb3VwLiBGb3IgdGhlIHRyZWF0bWVudCBncm91cCwgd2UgY2FuIGNhbGN1bGF0ZSB0aGUgaW5uZXIgdGVybSBmb3IgZWFjaCBmcmFtZSBzZXQgYW5kIHRoZW4gY29tYmluZSB0aGVpciBlZmZlY3RzLiBUaGUgdGVybSBpbiB0aGUgYnJhY2tldCBuZWVkcyBzcGVjaWFsIGF0dGVudGlvbiBhcyAkXHN1bV97Q157XHByaW1lfX1cbGVmdChwKENee1xwcmltZX0pcChZfFIsQ157XHByaW1lfSlccmlnaHQpJCBpcyBldmFsdWF0ZWQgZm9yIGJvdGggY2hvaWNlcyB3aGVuIGF0IGRpZmZlcmVudCAkUiQuIA0KDQpUaGUgcHJvY2VkdXJlIGZvciBldmFsdWF0aW9uIGlzIGxpa2UgdGhpcy4gRm9yICoqZWFjaCBmcmFtZSBzZXQgZ3JvdXAqKiwgd2UgdGhlIGZvbGxvd2luZzogDQoNCjEuIENhbGN1bGF0ZSB0aGlzIHRlcm0gJFxzdW1fe0Nee1xwcmltZX19XGxlZnQocChDXntccHJpbWV9KXAoWXxSLENee1xwcmltZX0pXHJpZ2h0KSQsIGF0ICRSPTEkIGFuZCAkUj0wJC4gVGhpcyBnaXZlcyB1c2UgJFAoWXxkbyhSKT0xKSQgYW5kICRQKFl8ZG8oUik9MCkkLg0KMi4gV2l0aCB0aGUgdmFsdWVzLCB3ZSBjYWxjdWxhdGUgdGhlIHZhbHVlIG9mICRQKFl8ZG8oQyk9MSkkIGFuZCAkUChZfGRvKEMpPTIpJCBhczoNCiQkXGJlZ2lue2VxbmFycmF5Kn0NClAoWXxkbyhDKT0xKSA9IFxzdW1fe1J9cChSfEM9MSlQKFl8ZG8oUl57XHByaW1lfSk9UilcXA0KUChZfGRvKEMpPTIpID0gXHN1bV97Un1wKFJ8Qz0yKVAoWXxkbyhSXntccHJpbWV9KT1SKQ0KXGVuZHtlcW5hcnJheSp9DQokJA0KDQpXaXRoIHZhbHVlcyBvZiAkUChZfGRvKEMpKSQsIHdlIGNvbWJpbmUgdGhlbSB0byBnZXQgdGhlIHRvdGFsIHBhdGggZWZmZWN0IGF0ICRGPWYkIChJIGlnbm9yZSB0aGUgJGRvKEYpPWYkIGFzICRGJCBpcyByYW5kb21pemVkKS4NCiQkUChZfGRvKEMpLEY9Zik9IFxzdW1fe2N9cChDPWN8Rj1mKVAoWXxkbyhDKT1jKSQkDQoNCkdpdmVuIHRoZSBwZXItZnJhbWUgc2V0IGdyb3VwIHZhbHVlcywgd2UgY29tYmluZSB0byBnZXQgdGhlIHRvdGFsIGluZGlyZWN0IGVmZmVjdCBmb3IgdGhlIHRyZWF0bWVudCBncm91cDogDQokJFAoWXxkbyhEKT0xLGRvKEMpKT1cc3VtX3tmfVAoRj1mfEQ9MSlQKFl8ZG8oQyksRj1mKSQk